From f3812079505bf98e9f51cef3bc55b7ebcc061a2e Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 14:23:22 +0200 Subject: [PATCH 1/9] moved code from all laskeMuutokset to same R file --- R/indMix.R | 204 ---------------------------------- R/laskeMuutokset12345.R | 236 ++++++++++++++++++++++++++++++++++++++++ 2 files changed, 236 insertions(+), 204 deletions(-) create mode 100644 R/laskeMuutokset12345.R 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 From f28bf5d7632ce9fa714673b17bfee37be56eb071 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 14:36:21 +0200 Subject: [PATCH 2/9] Translated laskeMuutokset1 --- R/laskeMuutokset12345.R | 71 ++++++++++++++++++++--------------------- R/laskeMuutokset4.R | 39 ---------------------- 2 files changed, 34 insertions(+), 76 deletions(-) delete mode 100644 R/laskeMuutokset4.R diff --git a/R/laskeMuutokset12345.R b/R/laskeMuutokset12345.R index 13c9196..410d6d0 100644 --- a/R/laskeMuutokset12345.R +++ b/R/laskeMuutokset12345.R @@ -40,51 +40,48 @@ laskeMuutokset4 <- function (osuus, osuusTaulu, omaFreqs, logml, 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. + # 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, ] - # npops = size(COUNTS,3); - # muutokset = LOGDIFF(ind,:); + i1 <- PARTITION[ind] + i1_logml <- POP_LOGML[i1] + muutokset[i1] <- 0 - # 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) - # 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 - # 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) - # 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) - # ni2 = length(i2); + COUNTS[, , i2] <- COUNTS[, , i2] + repmat(diffInCounts, c(1, 1, ni2)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] + repmat(diffInSumCounts, c(ni2, 1)) + new_i2_logml <- computePopulationLogml(i2, adjprior, priorTerm) + COUNTS[, , i2] <- COUNTS[, , i2] - repmat(diffInCounts, c(1, 1, ni2)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] - repmat(diffInSumCounts, c(ni2, 1)) - # 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; + muutokset[i2] <- new_i1_logml - i1_logml + new_i2_logml - i2_logml + LOGDIFF[ind, ] = muutokset return(list(muutokset = muutokset, diffInCounts = diffInCounts)) } diff --git a/R/laskeMuutokset4.R b/R/laskeMuutokset4.R deleted file mode 100644 index 56bcef3..0000000 --- a/R/laskeMuutokset4.R +++ /dev/null @@ -1,39 +0,0 @@ -#' @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) -} \ No newline at end of file From eba8705d9edd303d8b9f8707440f4f2ecd1390e5 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 15:32:34 +0200 Subject: [PATCH 3/9] Implemented basic functionality of setdiff --- R/setdiff.R | 13 +++++++++++++ tests/testthat/test-convertedBaseFunctions.R | 15 +++++++++++---- 2 files changed, 24 insertions(+), 4 deletions(-) create mode 100644 R/setdiff.R diff --git a/R/setdiff.R b/R/setdiff.R new file mode 100644 index 0000000..49e5633 --- /dev/null +++ b/R/setdiff.R @@ -0,0 +1,13 @@ +#' @title Set differences of two arrays +#' @description Loosely replicates the behavior of the homonym Matlab function +#' @param A first array +#' @param B second awway +#' @param legacy if `TRUE`, preserves the behavior of +#' @return +#' @author Waldir Leoncio +#' @export +setdiff <- function(A, B, legacy = FALSE) { + values <- sort(unique(A[is.na(match(A, B))])) + # browser() # TEMP + return(values) +} \ No newline at end of file diff --git a/tests/testthat/test-convertedBaseFunctions.R b/tests/testthat/test-convertedBaseFunctions.R index 11988f8..08e36cc 100644 --- a/tests/testthat/test-convertedBaseFunctions.R +++ b/tests/testthat/test-convertedBaseFunctions.R @@ -158,10 +158,10 @@ test_that("find works as expected", { }) test_that("sortrows works as expected", { - mx <- matrix(c(3, 2, 2, 1, 1, 10, 0, pi), 4) - expect_equal(sortrows(mx), matrix(c(1, 2, 2, 3, pi, 10, 0, 1), 4)) - expect_equal(sortrows(mx, 2), matrix(c(2, 3, 1, 2, 0, 1, pi, 10), 4)) - expect_equal(sortrows(mx, 1:2), mx[order(mx[, 1], mx[, 2]), ]) + mx <- matrix(c(3, 2, 2, 1, 1, 10, 0, pi), 4) + expect_equal(sortrows(mx), matrix(c(1, 2, 2, 3, pi, 10, 0, 1), 4)) + expect_equal(sortrows(mx, 2), matrix(c(2, 3, 1, 2, 0, 1, pi, 10), 4)) + expect_equal(sortrows(mx, 1:2), mx[order(mx[, 1], mx[, 2]), ]) }) test_that("cell works as expected", { @@ -217,4 +217,11 @@ test_that("nargin works correctly", { expect_equal(addme(13, 42), 55) expect_equal(addme(13), 26) expect_equal(addme(), 0) +}) + +test_that("setdiff works as expected", { + A <- c(3, 6, 2, 1, 5, 1, 1) + B <- c(2, 4, 6) + expect_equal(setdiff(A, B), list(c(1, 3, 5), 4, 1, 5)) + # TODO: add more examples from https://se.mathworks.com/help/matlab/ref/double.setdiff.html;jsessionid=0d8d42582d4d299b8224403899f1 }) \ No newline at end of file From 6b1a1910e3194bfc0919b815e9f4b95e40dac78e Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 15:57:15 +0200 Subject: [PATCH 4/9] Improvements to setdiff --- NAMESPACE | 1 + R/setdiff.R | 10 +++++++--- tests/testthat/test-convertedBaseFunctions.R | 19 ++++++++++++++++++- 3 files changed, 26 insertions(+), 4 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 2dd1c5e..51d5cf2 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -34,6 +34,7 @@ export(randdir) export(repmat) export(rivinSisaltamienMjonojenLkm) export(selvitaDigitFormat) +export(setdiff_MATLAB) export(simulateAllFreqs) export(simulateIndividuals) export(simuloiAlleeli) diff --git a/R/setdiff.R b/R/setdiff.R index 49e5633..af199a9 100644 --- a/R/setdiff.R +++ b/R/setdiff.R @@ -6,8 +6,12 @@ #' @return #' @author Waldir Leoncio #' @export -setdiff <- function(A, B, legacy = FALSE) { - values <- sort(unique(A[is.na(match(A, B))])) - # browser() # TEMP +setdiff_MATLAB <- function(A, B, legacy = FALSE) { + if (is(A, "numeric") & is(B, "numeric")) { + values <- sort(unique(A[is.na(match(A, B))])) + } else if (is(A, "data.frame") & is(B, "data.frame")) { + stop("Not implemented for data frames") + } + # TODO: add support for indices (if necessary) return(values) } \ No newline at end of file diff --git a/tests/testthat/test-convertedBaseFunctions.R b/tests/testthat/test-convertedBaseFunctions.R index 08e36cc..01f6e2e 100644 --- a/tests/testthat/test-convertedBaseFunctions.R +++ b/tests/testthat/test-convertedBaseFunctions.R @@ -222,6 +222,23 @@ test_that("nargin works correctly", { test_that("setdiff works as expected", { A <- c(3, 6, 2, 1, 5, 1, 1) B <- c(2, 4, 6) - expect_equal(setdiff(A, B), list(c(1, 3, 5), 4, 1, 5)) + C <- c(1, 3, 5) + expect_equal(setdiff(A, B), C) + A <- data.frame( + Var1 = 1:5, + Var2 = LETTERS[1:5], + Var3 = c(FALSE, TRUE, FALSE, TRUE, FALSE) + ) + B <- data.frame( + Var1 = seq(1, 9, by = 2), + Var2 = LETTERS[seq(1, 9, by = 2)], + Var3 = rep(FALSE, 5) + ) + C <- data.frame( + Var1 = c(2, 4), + Var2 = c('B', 'D'), + Var3 = c(TRUE, TRUE) + ) + expect_equal(setdiff(A, B), C) # TODO: add more examples from https://se.mathworks.com/help/matlab/ref/double.setdiff.html;jsessionid=0d8d42582d4d299b8224403899f1 }) \ No newline at end of file From 637df694bb2559e1346160a623b07f20412475c1 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 15:57:30 +0200 Subject: [PATCH 5/9] Updated docs --- NAMESPACE | 1 + man/admixture_initialization.Rd | 11 +++++++++++ man/laskeMuutokset4.Rd | 4 ++-- man/max_MATLAB.Rd | 22 ++++++++++++++++++++++ man/min_MATLAB.Rd | 18 +++++++++++++++++- man/setdiff_MATLAB.Rd | 24 ++++++++++++++++++++++++ 6 files changed, 77 insertions(+), 3 deletions(-) create mode 100644 man/admixture_initialization.Rd create mode 100644 man/max_MATLAB.Rd create mode 100644 man/setdiff_MATLAB.Rd diff --git a/NAMESPACE b/NAMESPACE index 51d5cf2..94690c0 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -23,6 +23,7 @@ export(linkage) export(logml2String) export(lueGenePopData) export(lueNimi) +export(max_MATLAB) export(min_MATLAB) export(noIndex) export(ownNum2Str) diff --git a/man/admixture_initialization.Rd b/man/admixture_initialization.Rd new file mode 100644 index 0000000..ae2c9bc --- /dev/null +++ b/man/admixture_initialization.Rd @@ -0,0 +1,11 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/admixture_initialization.R +\name{admixture_initialization} +\alias{admixture_initialization} +\title{Seuraavat kolme funktiota liittyvat alkupartition muodostamiseen.} +\usage{ +admixture_initialization(data_matrix, nclusters, Z) +} +\description{ +Seuraavat kolme funktiota liittyvat alkupartition muodostamiseen. +} diff --git a/man/laskeMuutokset4.Rd b/man/laskeMuutokset4.Rd index 7c16a19..5c41852 100644 --- a/man/laskeMuutokset4.Rd +++ b/man/laskeMuutokset4.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/laskeMuutokset4.R +% Please edit documentation in R/laskeMuutokset12345.R \name{laskeMuutokset4} \alias{laskeMuutokset4} \title{Calculate changes?} @@ -20,6 +20,6 @@ laskeMuutokset4(osuus, osuusTaulu, omaFreqs, logml, COUNTS = matrix(0)) \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 +todennäköisyysmassaa populaatioon j. Mikäli populaatiossa i ei ole mitään siirrettävää, on vastaavassa kohdassa rivi nollia. } diff --git a/man/max_MATLAB.Rd b/man/max_MATLAB.Rd new file mode 100644 index 0000000..d7fbdf7 --- /dev/null +++ b/man/max_MATLAB.Rd @@ -0,0 +1,22 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/min_max_MATLAB.R +\name{max_MATLAB} +\alias{max_MATLAB} +\title{Maximum (MATLAB version)} +\usage{ +max_MATLAB(X, indices = TRUE) +} +\arguments{ +\item{X}{matrix} + +\item{indices}{return indices?} +} +\value{ +Either a list or a vector +} +\description{ +Finds the minimum value for each column of a matrix, potentially returning the indices instead +} +\author{ +Waldir Leoncio +} diff --git a/man/min_MATLAB.Rd b/man/min_MATLAB.Rd index bd1113e..9bb1166 100644 --- a/man/min_MATLAB.Rd +++ b/man/min_MATLAB.Rd @@ -1,11 +1,15 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/min.R, R/min_MATLAB.R +% Please edit documentation in R/min.R, R/min_MATLAB.R, R/min_max_MATLAB.R \name{min_MATLAB} \alias{min_MATLAB} \title{Minimum (MATLAB version)} \usage{ min_MATLAB(X, indices = TRUE) +min_MATLAB(X, indices = TRUE) + +min_MATLAB(X, indices = TRUE) + min_MATLAB(X, indices = TRUE) } \arguments{ @@ -16,15 +20,27 @@ min_MATLAB(X, indices = TRUE) \value{ Either a list or a vector +Either a list or a vector + +Either a list or a vector + Either a list or a vector } \description{ Finds the minimum value for each column of a matrix, potentially returning the indices instead +Finds the minimum value for each column of a matrix, potentially returning the indices instead + +Finds the minimum value for each column of a matrix, potentially returning the indices instead + Finds the minimum value for each column of a matrix, potentially returning the indices instead } \author{ Waldir Leoncio +Waldir Leoncio + +Waldir Leoncio + Waldir Leoncio } diff --git a/man/setdiff_MATLAB.Rd b/man/setdiff_MATLAB.Rd new file mode 100644 index 0000000..681e9f9 --- /dev/null +++ b/man/setdiff_MATLAB.Rd @@ -0,0 +1,24 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/setdiff.R +\name{setdiff_MATLAB} +\alias{setdiff_MATLAB} +\title{Set differences of two arrays} +\usage{ +setdiff_MATLAB(A, B, legacy = FALSE) +} +\arguments{ +\item{A}{first array} + +\item{B}{second awway} + +\item{legacy}{if `TRUE`, preserves the behavior of} +} +\value{ + +} +\description{ +Loosely replicates the behavior of the homonym Matlab function +} +\author{ +Waldir Leoncio +} From 74e11205b72092508a13a18125f1bdf5a9704fcd Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 15:58:08 +0200 Subject: [PATCH 6/9] Fixed syntax --- R/addToSummary.R | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/R/addToSummary.R b/R/addToSummary.R index 1c79147..4affbe9 100644 --- a/R/addToSummary.R +++ b/R/addToSummary.R @@ -1,4 +1,4 @@ -addToSummary <- funciton(logml, partitionSummary, worstIndex) { +addToSummary <- function(logml, partitionSummary, worstIndex) { # Tiedet��n, ett� annettu logml on isompi kuin huonoin arvo # partitionSummary taulukossa. Jos partitionSummary:ss� ei viel� ole # annettua logml arvoa, niin lis�t��n worstIndex:in kohtaan uusi logml ja @@ -6,13 +6,13 @@ addToSummary <- funciton(logml, partitionSummary, worstIndex) { apu <- find(abs(partitionSummary[, 2] - logml) < 1e-5) if (isempty(apu)) { - # Nyt l�ydetty partitio ei ole viel� kirjattuna summaryyn. - npops <- length(unique(PARTITION)) - partitionSummary[worstIndex, 1] <- npops - partitionSummary[worstIndex, 2] <- logml - added <- 1 + # Nyt l�ydetty partitio ei ole viel� kirjattuna summaryyn. + npops <- length(unique(PARTITION)) + partitionSummary[worstIndex, 1] <- npops + partitionSummary[worstIndex, 2] <- logml + added <- 1 } else { - added <- 0 + added <- 0 } return(list(partitionSummary = partitionSummary, added = added)) } \ No newline at end of file From 448373796cbf80d8498a2314e6f42552fa4da001 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 15:58:17 +0200 Subject: [PATCH 7/9] Added FIXME --- tests/testthat/test-convertedBaseFunctions.R | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/testthat/test-convertedBaseFunctions.R b/tests/testthat/test-convertedBaseFunctions.R index 01f6e2e..5623e38 100644 --- a/tests/testthat/test-convertedBaseFunctions.R +++ b/tests/testthat/test-convertedBaseFunctions.R @@ -164,6 +164,7 @@ test_that("sortrows works as expected", { expect_equal(sortrows(mx, 1:2), mx[order(mx[, 1], mx[, 2]), ]) }) +# FIXME: failing tests test_that("cell works as expected", { expect_equal(cell(0), array(dim = c(0, 0))) expect_equal(cell(1), array(dim = c(1, 1))) From 9a38c59ab98c0224aebd573e0413d4456441857b Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 16:12:22 +0200 Subject: [PATCH 8/9] Translated laskeMuutokset2 and 3 --- R/laskeMuutokset12345.R | 149 +++++++++++++++++++--------------------- 1 file changed, 71 insertions(+), 78 deletions(-) diff --git a/R/laskeMuutokset12345.R b/R/laskeMuutokset12345.R index 410d6d0..3a4f3ae 100644 --- a/R/laskeMuutokset12345.R +++ b/R/laskeMuutokset12345.R @@ -70,7 +70,7 @@ laskeMuutokset <- function(ind, globalRows, data, adjprior, priorTerm) { i2 <- find(muutokset == -Inf) # Etsit��n populaatiot jotka muuttuneet viime kerran j�lkeen. i2 <- setdiff(i2, i1) - i2_logml <- POP_LOGML(i2) + i2_logml <- POP_LOGML[i2] ni2 <- length(i2) @@ -85,107 +85,100 @@ laskeMuutokset <- function(ind, globalRows, data, adjprior, priorTerm) { 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); + npops <- size(COUNTS, 3) + muutokset <- zeros(npops, 1) - # i1_logml = POP_LOGML(i1); + i1_logml <- POP_LOGML[i1] - # inds = find(PARTITION==i1); - # ninds = length(inds); + inds <- find(PARTITION == i1) + ninds <- length(inds) - # if ninds==0 - # diffInCounts = zeros(size(COUNTS,1), size(COUNTS,2)); - # return; - # end + if (ninds == 0) { + diffInCounts <- zeros(size(COUNTS, 1), size(COUNTS, 2)) + return() + } - # rows = []; - # for i = 1:ninds - # ind = inds(i); - # lisa = globalRows(ind,1):globalRows(ind,2); - # rows = [rows; lisa']; - # %rows = [rows; globalRows{ind}']; - # end + rows = list() + for (i in 1:ninds) { + ind <- inds(i) + lisa <- globalRows(ind, 1):globalRows(ind, 2) + rows <- c(rows, t(lisa)) + } - # diffInCounts = computeDiffInCounts(rows', size(COUNTS,1), size(COUNTS,2), data); - # diffInSumCounts = sum(diffInCounts); + diffInCounts <- computeDiffInCounts( + t(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; + 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); + i2 <- c(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]); + COUNTS[, , i2] <- COUNTS[, , i2] + repmat(diffInCounts, c(1, 1, npops - 1)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] + repmat(diffInSumCounts, c(npops - 1, 1)) + new_i2_logml <- computePopulationLogml(i2, adjprior, priorTerm) + COUNTS[, , i2] <- COUNTS[, , i2] - repmat(diffInCounts, c(1, 1, npops - 1)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] - repmat(diffInSumCounts, c(npops - 1, 1)) - # muutokset(i2) = new_i1_logml - i1_logml ... - # + new_i2_logml - i2_logml; - return(list(muutokset = muutokset, diffInCounts = diffInCounts)) + 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. + # 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); + 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); + i1_logml = POP_LOGML[i1] + for (pop2 in 1:npops2) { + inds <- inds2[find(T2==pop2)] + ninds <- length(inds); + if (ninds > 0) { + rows <- list() + for (i in 1:ninds) { + ind <- inds[i] + lisa <- globalRows[ind, 1]:globalRows[ind, 2] + rows <- c(rows, t(lisa)) + } + diffInCounts <- computeDiffInCounts( + t(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; + 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)'; + i2 <- c(1:i1-1, i1+1:npops) + i2_logml <- t(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]); + COUNTS[, , i2] <- COUNTS[, , i2] + repmat(diffInCounts, c(1, 1, npops - 1)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] + repmat(diffInSumCounts, c(npops - 1, 1)) + new_i2_logml <- t(computePopulationLogml(i2, adjprior, priorTerm)) + COUNTS[, , i2] <- COUNTS[, , i2] - repmat(diffInCounts, c(1, 1, npops - 1)) + SUMCOUNTS[i2, ] <- SUMCOUNTS[i2, ] - repmat(diffInSumCounts, c(npops - 1, 1)) - # muutokset(pop2,i2) = new_i1_logml - i1_logml ... - # + new_i2_logml - i2_logml; - # end - # end - return(muutokset) + muutokset[pop2, i2] <- new_i1_logml - i1_logml + new_i2_logml - i2_logml + } + } + return(muutokset) } laskeMuutokset5 <- function(inds, globalRows, data, adjprior, priorTerm, i1, i2) { From 6d6d8c557dc5eebc3a60019e8f23a11c00ef5460 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Mon, 19 Oct 2020 16:22:53 +0200 Subject: [PATCH 9/9] Added lakseMuutokset5 --- R/laskeMuutokset12345.R | 67 +++++++++++++++++++++-------------------- 1 file changed, 34 insertions(+), 33 deletions(-) diff --git a/R/laskeMuutokset12345.R b/R/laskeMuutokset12345.R index 3a4f3ae..b642a98 100644 --- a/R/laskeMuutokset12345.R +++ b/R/laskeMuutokset12345.R @@ -182,45 +182,46 @@ laskeMuutokset3 <- function(T2, inds2, globalRows, data, adjprior, priorTerm, i1 } 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�. + # 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) - # ninds = length(inds); - # muutokset = zeros(ninds,1); + i1_logml <- POP_LOGML[i1] + i2_logml <- POP_LOGML[i2] - # i1_logml = POP_LOGML(i1); - # i2_logml = POP_LOGML(i2); + for (i in 1:ninds) { + ind <- inds[i] + if (PARTITION[ind] == i1) { + pop1 <- i1 #mist� + pop2 <- i2 #mihin + } else { + pop1 <- i2 + pop2 <- i1 + } + rows <- globalRows[ind, 1]:globalRows[ind, 2] + diffInCounts <- computeDiffInCounts( + rows, size(COUNTS, 1), size(COUNTS, 2), data + ) + diffInSumCounts <- sum(diffInCounts) - # 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 - # COUNTS(:,:,pop1) = COUNTS(:,:,pop1)+diffInCounts; - # SUMCOUNTS(pop1,:) = SUMCOUNTS(pop1,:)+diffInSumCounts; - # COUNTS(:,:,pop2) = COUNTS(:,:,pop2)-diffInCounts; - # SUMCOUNTS(pop2,:) = SUMCOUNTS(pop2,:)-diffInSumCounts; - # end + new_logmls <- computePopulationLogml(c(i1, i2), adjprior, priorTerm) + muutokset[i] <- sum(new_logmls) - # muutokset = muutokset - i1_logml - i2_logml; - return(muutokset) + COUNTS[, , pop1] <- COUNTS[, , pop1] + diffInCounts + SUMCOUNTS[pop1, ] <- SUMCOUNTS[pop1, ] + diffInSumCounts + COUNTS[, , pop2] <- COUNTS[, , pop2] - diffInCounts + SUMCOUNTS[pop2, ] <- SUMCOUNTS[pop2, ] - diffInSumCounts + } + + muutokset <- muutokset - i1_logml - i2_logml + return(muutokset) } \ No newline at end of file