diff --git a/R/indMix.R b/R/indMix.R index 62d488a..91fe4dd 100644 --- a/R/indMix.R +++ b/R/indMix.R @@ -678,57 +678,6 @@ indMix <- function(c, npops, dispText) { # %------------------------------------------------------------------------------------ -# function [muutokset, diffInCounts] = ... -# laskeMuutokset(ind, globalRows, data, adjprior, priorTerm) -# % Palauttaa npops*1 taulun, jossa i:s alkio kertoo, mik� olisi -# % muutos logml:ss�, mik�li yksil� ind siirret��n koriin i. -# % diffInCounts on poistettava COUNTS:in siivusta i1 ja lis�tt�v� -# % COUNTS:in siivuun i2, mik�li muutos toteutetaan. -# % -# % Lis�ys 25.9.2007: -# % Otettu k�ytt��n globaali muuttuja LOGDIFF, johon on tallennettu muutokset -# % logml:ss� siirrett�ess� yksil�it� toisiin populaatioihin. - -# global COUNTS; global SUMCOUNTS; -# global PARTITION; global POP_LOGML; -# global LOGDIFF; - -# npops = size(COUNTS,3); -# muutokset = LOGDIFF(ind,:); - -# i1 = PARTITION(ind); -# i1_logml = POP_LOGML(i1); -# muutokset(i1) = 0; - -# rows = globalRows(ind,1):globalRows(ind,2); -# diffInCounts = computeDiffInCounts(rows, size(COUNTS,1), size(COUNTS,2), data); -# diffInSumCounts = sum(diffInCounts); - -# COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; -# new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); -# COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; - -# i2 = find(muutokset==-Inf); % Etsit��n populaatiot jotka muuttuneet viime kerran j�lkeen. -# i2 = setdiff(i2,i1); -# i2_logml = POP_LOGML(i2); - -# ni2 = length(i2); - -# COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 ni2]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[ni2 1]); -# new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm); -# COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 ni2]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[ni2 1]); - -# muutokset(i2) = new_i1_logml - i1_logml ... -# + new_i2_logml - i2_logml; -# LOGDIFF(ind,:) = muutokset; - - -# %---------------------------------------------------------------------- - # function diffInCounts = computeDiffInCounts(rows, max_noalle, nloci, data) # % Muodostaa max_noalle*nloci taulukon, jossa on niiden alleelien @@ -781,62 +730,6 @@ indMix <- function(c, npops, dispText) { # %-------------------------------------------------------------------------- # %-- -# %------------------------------------------------------------------------------------ - - -# function [muutokset, diffInCounts] = laskeMuutokset2( ... -# i1, globalRows, data, adjprior, priorTerm); -# % Palauttaa npops*1 taulun, jossa i:s alkio kertoo, mik� olisi -# % muutos logml:ss�, mik�li korin i1 kaikki yksil�t siirret��n -# % koriin i. - -# global COUNTS; global SUMCOUNTS; -# global PARTITION; global POP_LOGML; -# npops = size(COUNTS,3); -# muutokset = zeros(npops,1); - -# i1_logml = POP_LOGML(i1); - -# inds = find(PARTITION==i1); -# ninds = length(inds); - -# if ninds==0 -# diffInCounts = zeros(size(COUNTS,1), size(COUNTS,2)); -# return; -# end - -# rows = []; -# for i = 1:ninds -# ind = inds(i); -# lisa = globalRows(ind,1):globalRows(ind,2); -# rows = [rows; lisa']; -# %rows = [rows; globalRows{ind}']; -# end - -# diffInCounts = computeDiffInCounts(rows', size(COUNTS,1), size(COUNTS,2), data); -# diffInSumCounts = sum(diffInCounts); - -# COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; -# new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); -# COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; - -# i2 = [1:i1-1 , i1+1:npops]; -# i2_logml = POP_LOGML(i2); - -# COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 npops-1]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[npops-1 1]); -# new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm); -# COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 npops-1]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[npops-1 1]); - -# muutokset(i2) = new_i1_logml - i1_logml ... -# + new_i2_logml - i2_logml; - - -# %--------------------------------------------------------------------------------- - # function updateGlobalVariables2( ... # i1, i2, diffInCounts, adjprior, priorTerm); @@ -864,103 +757,6 @@ indMix <- function(c, npops, dispText) { # inx = [find(PARTITION==i1); find(PARTITION==i2)]; # LOGDIFF(inx,:) = -Inf; - -# %-------------------------------------------------------------------------- -# %---- - -# function muutokset = laskeMuutokset3(T2, inds2, globalRows, ... -# data, adjprior, priorTerm, i1) -# % Palauttaa length(unique(T2))*npops taulun, jossa (i,j):s alkio -# % kertoo, mik� olisi muutos logml:ss�, jos populaation i1 osapopulaatio -# % inds2(find(T2==i)) siirret��n koriin j. - -# global COUNTS; global SUMCOUNTS; -# global PARTITION; global POP_LOGML; -# npops = size(COUNTS,3); -# npops2 = length(unique(T2)); -# muutokset = zeros(npops2, npops); - -# i1_logml = POP_LOGML(i1); -# for pop2 = 1:npops2 -# inds = inds2(find(T2==pop2)); -# ninds = length(inds); -# if ninds>0 -# rows = []; -# for i = 1:ninds -# ind = inds(i); -# lisa = globalRows(ind,1):globalRows(ind,2); -# rows = [rows; lisa']; -# %rows = [rows; globalRows{ind}']; -# end -# diffInCounts = computeDiffInCounts(rows', size(COUNTS,1), size(COUNTS,2), data); -# diffInSumCounts = sum(diffInCounts); - -# COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; -# new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); -# COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; -# SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; - -# i2 = [1:i1-1 , i1+1:npops]; -# i2_logml = POP_LOGML(i2)'; - -# COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 npops-1]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[npops-1 1]); -# new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm)'; -# COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 npops-1]); -# SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[npops-1 1]); - -# muutokset(pop2,i2) = new_i1_logml - i1_logml ... -# + new_i2_logml - i2_logml; -# end -# end - -# %------------------------------------------------------------------------------------ - -# function muutokset = laskeMuutokset5(inds, globalRows, data, adjprior, ... -# priorTerm, i1, i2) - -# % Palauttaa length(inds)*1 taulun, jossa i:s alkio kertoo, mik� olisi -# % muutos logml:ss�, mik�li yksil� i vaihtaisi koria i1:n ja i2:n v�lill�. - -# global COUNTS; global SUMCOUNTS; -# global PARTITION; global POP_LOGML; - -# ninds = length(inds); -# muutokset = zeros(ninds,1); - -# i1_logml = POP_LOGML(i1); -# i2_logml = POP_LOGML(i2); - -# for i = 1:ninds -# ind = inds(i); -# if PARTITION(ind)==i1 -# pop1 = i1; %mist� -# pop2 = i2; %mihin -# else -# pop1 = i2; -# pop2 = i1; -# end -# rows = globalRows(ind,1):globalRows(ind,2); -# diffInCounts = computeDiffInCounts(rows, size(COUNTS,1), size(COUNTS,2), data); -# diffInSumCounts = sum(diffInCounts); - -# COUNTS(:,:,pop1) = COUNTS(:,:,pop1)-diffInCounts; -# SUMCOUNTS(pop1,:) = SUMCOUNTS(pop1,:)-diffInSumCounts; -# COUNTS(:,:,pop2) = COUNTS(:,:,pop2)+diffInCounts; -# SUMCOUNTS(pop2,:) = SUMCOUNTS(pop2,:)+diffInSumCounts; - -# new_logmls = computePopulationLogml([i1 i2], adjprior, priorTerm); -# muutokset(i) = sum(new_logmls); - -# COUNTS(:,:,pop1) = COUNTS(:,:,pop1)+diffInCounts; -# SUMCOUNTS(pop1,:) = SUMCOUNTS(pop1,:)+diffInSumCounts; -# COUNTS(:,:,pop2) = COUNTS(:,:,pop2)-diffInCounts; -# SUMCOUNTS(pop2,:) = SUMCOUNTS(pop2,:)-diffInSumCounts; -# end - -# muutokset = muutokset - i1_logml - i2_logml; - # %------------------------------------------------------------------------------------ diff --git a/R/laskeMuutokset12345.R b/R/laskeMuutokset12345.R new file mode 100644 index 0000000..13c9196 --- /dev/null +++ b/R/laskeMuutokset12345.R @@ -0,0 +1,236 @@ +#' @title Calculate changes? +#' @description Palauttaa npops*npops taulun, jonka alkio (i,j) kertoo, mik?on +#' muutos logml:ss? mikäli populaatiosta i siirretään osuuden verran +#' todennäköisyysmassaa populaatioon j. Mikäli populaatiossa i ei ole mitään +#' siirrettävää, on vastaavassa kohdassa rivi nollia. +#' @param osuus Percentages? +#' @param omaFreqs own Freqs? +#' @param osuusTaulu Percentage table? +#' @param logml log maximum likelihood +#' @param COUNTS COUNTS +#' @export +laskeMuutokset4 <- function (osuus, osuusTaulu, omaFreqs, logml, + COUNTS = matrix(0)) { + npops <- ifelse(is.na(dim(COUNTS)[3]), 1, dim(COUNTS)[3]) + notEmpty <- which(osuusTaulu > 0.005) + muutokset <- zeros(npops) + empties <- !notEmpty + + for (i1 in notEmpty) { + osuusTaulu[i1] <- osuusTaulu[i1] - osuus + for (i2 in c(colon(1, i1 - 1), colon(i1 + 1, npops))) { + osuusTaulu[i2] <- osuusTaulu[i2] + osuus + loggis <- computeIndLogml(omaFreqs, osuusTaulu) + + # Work around Matlab OOB bug + if (i1 > nrow(muutokset)) { + muutokset <- rbind(muutokset, muutokset * 0) + } + if (i2 > ncol(muutokset)) { + muutokset <- cbind(muutokset, muutokset * 0) + } + + muutokset[i1, i2] <- loggis - logml + osuusTaulu[i2] <- osuusTaulu[i2] - osuus + } + osuusTaulu[i1] <- osuusTaulu[i1] + osuus + } + return (muutokset) +} + + +laskeMuutokset <- function(ind, globalRows, data, adjprior, priorTerm) { + # % Palauttaa npops*1 taulun, jossa i:s alkio kertoo, mik� olisi + # % muutos logml:ss�, mik�li yksil� ind siirret��n koriin i. + # % diffInCounts on poistettava COUNTS:in siivusta i1 ja lis�tt�v� + # % COUNTS:in siivuun i2, mik�li muutos toteutetaan. + # % + # % Lis�ys 25.9.2007: + # % Otettu k�ytt��n globaali muuttuja LOGDIFF, johon on tallennettu muutokset + # % logml:ss� siirrett�ess� yksil�it� toisiin populaatioihin. + + # global COUNTS; global SUMCOUNTS; + # global PARTITION; global POP_LOGML; + # global LOGDIFF; + + # npops = size(COUNTS,3); + # muutokset = LOGDIFF(ind,:); + + # i1 = PARTITION(ind); + # i1_logml = POP_LOGML(i1); + # muutokset(i1) = 0; + + # rows = globalRows(ind,1):globalRows(ind,2); + # diffInCounts = computeDiffInCounts(rows, size(COUNTS,1), size(COUNTS,2), data); + # diffInSumCounts = sum(diffInCounts); + + # COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; + # new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); + # COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; + + # i2 = find(muutokset==-Inf); % Etsit��n populaatiot jotka muuttuneet viime kerran j�lkeen. + # i2 = setdiff(i2,i1); + # i2_logml = POP_LOGML(i2); + + # ni2 = length(i2); + + # COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 ni2]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[ni2 1]); + # new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm); + # COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 ni2]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[ni2 1]); + + # muutokset(i2) = new_i1_logml - i1_logml ... + # + new_i2_logml - i2_logml; + # LOGDIFF(ind,:) = muutokset; + return(list(muutokset = muutokset, diffInCounts = diffInCounts)) +} + +# %------------------------------------------------------------------------------------ + + +laskeMuutokset2 <- function(i1, globalRows, data, adjprior, priorTerm) { + # % Palauttaa npops*1 taulun, jossa i:s alkio kertoo, mik� olisi + # % muutos logml:ss�, mik�li korin i1 kaikki yksil�t siirret��n + # % koriin i. + + # global COUNTS; global SUMCOUNTS; + # global PARTITION; global POP_LOGML; + # npops = size(COUNTS,3); + # muutokset = zeros(npops,1); + + # i1_logml = POP_LOGML(i1); + + # inds = find(PARTITION==i1); + # ninds = length(inds); + + # if ninds==0 + # diffInCounts = zeros(size(COUNTS,1), size(COUNTS,2)); + # return; + # end + + # rows = []; + # for i = 1:ninds + # ind = inds(i); + # lisa = globalRows(ind,1):globalRows(ind,2); + # rows = [rows; lisa']; + # %rows = [rows; globalRows{ind}']; + # end + + # diffInCounts = computeDiffInCounts(rows', size(COUNTS,1), size(COUNTS,2), data); + # diffInSumCounts = sum(diffInCounts); + + # COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; + # new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); + # COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; + + # i2 = [1:i1-1 , i1+1:npops]; + # i2_logml = POP_LOGML(i2); + + # COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 npops-1]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[npops-1 1]); + # new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm); + # COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 npops-1]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[npops-1 1]); + + # muutokset(i2) = new_i1_logml - i1_logml ... + # + new_i2_logml - i2_logml; + return(list(muutokset = muutokset, diffInCounts = diffInCounts)) +} + + +laskeMuutokset3 <- function(T2, inds2, globalRows, data, adjprior, priorTerm, i1) { + # % Palauttaa length(unique(T2))*npops taulun, jossa (i,j):s alkio + # % kertoo, mik� olisi muutos logml:ss�, jos populaation i1 osapopulaatio + # % inds2(find(T2==i)) siirret��n koriin j. + + # global COUNTS; global SUMCOUNTS; + # global PARTITION; global POP_LOGML; + # npops = size(COUNTS,3); + # npops2 = length(unique(T2)); + # muutokset = zeros(npops2, npops); + + # i1_logml = POP_LOGML(i1); + # for pop2 = 1:npops2 + # inds = inds2(find(T2==pop2)); + # ninds = length(inds); + # if ninds>0 + # rows = []; + # for i = 1:ninds + # ind = inds(i); + # lisa = globalRows(ind,1):globalRows(ind,2); + # rows = [rows; lisa']; + # %rows = [rows; globalRows{ind}']; + # end + # diffInCounts = computeDiffInCounts(rows', size(COUNTS,1), size(COUNTS,2), data); + # diffInSumCounts = sum(diffInCounts); + + # COUNTS(:,:,i1) = COUNTS(:,:,i1)-diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)-diffInSumCounts; + # new_i1_logml = computePopulationLogml(i1, adjprior, priorTerm); + # COUNTS(:,:,i1) = COUNTS(:,:,i1)+diffInCounts; + # SUMCOUNTS(i1,:) = SUMCOUNTS(i1,:)+diffInSumCounts; + + # i2 = [1:i1-1 , i1+1:npops]; + # i2_logml = POP_LOGML(i2)'; + + # COUNTS(:,:,i2) = COUNTS(:,:,i2)+repmat(diffInCounts, [1 1 npops-1]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)+repmat(diffInSumCounts,[npops-1 1]); + # new_i2_logml = computePopulationLogml(i2, adjprior, priorTerm)'; + # COUNTS(:,:,i2) = COUNTS(:,:,i2)-repmat(diffInCounts, [1 1 npops-1]); + # SUMCOUNTS(i2,:) = SUMCOUNTS(i2,:)-repmat(diffInSumCounts,[npops-1 1]); + + # muutokset(pop2,i2) = new_i1_logml - i1_logml ... + # + new_i2_logml - i2_logml; + # end + # end + return(muutokset) +} + +laskeMuutokset5 <- function(inds, globalRows, data, adjprior, priorTerm, i1, i2) { + # % Palauttaa length(inds)*1 taulun, jossa i:s alkio kertoo, mik� olisi + # % muutos logml:ss�, mik�li yksil� i vaihtaisi koria i1:n ja i2:n v�lill�. + + # global COUNTS; global SUMCOUNTS; + # global PARTITION; global POP_LOGML; + + # ninds = length(inds); + # muutokset = zeros(ninds,1); + + # i1_logml = POP_LOGML(i1); + # i2_logml = POP_LOGML(i2); + + # for i = 1:ninds + # ind = inds(i); + # if PARTITION(ind)==i1 + # pop1 = i1; %mist� + # pop2 = i2; %mihin + # else + # pop1 = i2; + # pop2 = i1; + # end + # rows = globalRows(ind,1):globalRows(ind,2); + # diffInCounts = computeDiffInCounts(rows, size(COUNTS,1), size(COUNTS,2), data); + # diffInSumCounts = sum(diffInCounts); + + # COUNTS(:,:,pop1) = COUNTS(:,:,pop1)-diffInCounts; + # SUMCOUNTS(pop1,:) = SUMCOUNTS(pop1,:)-diffInSumCounts; + # COUNTS(:,:,pop2) = COUNTS(:,:,pop2)+diffInCounts; + # SUMCOUNTS(pop2,:) = SUMCOUNTS(pop2,:)+diffInSumCounts; + + # new_logmls = computePopulationLogml([i1 i2], adjprior, priorTerm); + # muutokset(i) = sum(new_logmls); + + # COUNTS(:,:,pop1) = COUNTS(:,:,pop1)+diffInCounts; + # SUMCOUNTS(pop1,:) = SUMCOUNTS(pop1,:)+diffInSumCounts; + # COUNTS(:,:,pop2) = COUNTS(:,:,pop2)-diffInCounts; + # SUMCOUNTS(pop2,:) = SUMCOUNTS(pop2,:)-diffInSumCounts; + # end + + # muutokset = muutokset - i1_logml - i2_logml; + return(muutokset) +} \ No newline at end of file