Code cleanup

This commit is contained in:
Waldir Leoncio 2020-07-14 11:18:20 +02:00
parent c2bbe701be
commit a78a9c8a7d

View file

@ -366,22 +366,6 @@ greedyMix <- function(
# global PARTITION; PARTITION = []; # global PARTITION; PARTITION = [];
# global POP_LOGML; POP_LOGML = []; # global POP_LOGML; POP_LOGML = [];
# %-------------------------------------------------------------------------------------
# function rows = computeRows(rowsFromInd, inds, ninds)
# % On annettu yksil<69>t inds. Funktio palauttaa vektorin, joka
# % sis<69>lt<6C><74> niiden rivien numerot, jotka sis<69>lt<6C>v<EFBFBD>t yksil<69>iden
# % dataa.
# rows = inds(:, ones(1,rowsFromInd));
# rows = rows*rowsFromInd;
# miinus = repmat(rowsFromInd-1 : -1 : 0, [ninds 1]);
# rows = rows - miinus;
# rows = reshape(rows', [1,rowsFromInd*ninds]);
# %-------------------------------------------------------------------------- # %--------------------------------------------------------------------------
@ -1210,270 +1194,6 @@ greedyMix <- function(
# end # end
# end # end
# %-------------------------------------------------------------------
# function changesInLogml = writeMixtureInfo(logml, rowsFromInd, data, adjprior, ...
# priorTerm, outPutFile, inputFile, partitionSummary, popnames, fixedK)
# global PARTITION;
# global COUNTS;
# global SUMCOUNTS;
# global LOGDIFF;
# changesInLogml = [];
# ninds = size(data,1)/rowsFromInd;
# npops = size(COUNTS,3);
# names = (size(popnames,1) == ninds); %Tarkistetaan ett?nimet viittaavat yksil<69>ihin
# if length(outPutFile)>0
# fid = fopen(outPutFile,'a');
# else
# fid = -1;
# diary('baps4_output.baps'); % save in text anyway.
# end
# dispLine;
# disp('RESULTS OF INDIVIDUAL LEVEL MIXTURE ANALYSIS:');
# disp(['Data file: ' inputFile]);
# disp(['Model: independent']);
# disp(['Number of clustered individuals: ' ownNum2Str(ninds)]);
# disp(['Number of groups in optimal partition: ' ownNum2Str(npops)]);
# disp(['Log(marginal likelihood) of optimal partition: ' ownNum2Str(logml)]);
# disp(' ');
# if (fid ~= -1)
# fprintf(fid,'%s \n', ['RESULTS OF INDIVIDUAL LEVEL MIXTURE ANALYSIS:']); fprintf(fid,'\n');
# fprintf(fid,'%s \n', ['Data file: ' inputFile]); fprintf(fid,'\n');
# fprintf(fid,'%s \n', ['Number of clustered individuals: ' ownNum2Str(ninds)]); fprintf(fid,'\n');
# fprintf(fid,'%s \n', ['Number of groups in optimal partition: ' ownNum2Str(npops)]); fprintf(fid,'\n');
# fprintf(fid,'%s \n', ['Log(marginal likelihood) of optimal partition: ' ownNum2Str(logml)]); fprintf(fid,'\n');
# end
# cluster_count = length(unique(PARTITION));
# disp(['Best Partition: ']);
# if (fid ~= -1)
# fprintf(fid,'%s \n',['Best Partition: ']); fprintf(fid,'\n');
# end
# for m=1:cluster_count
# indsInM = find(PARTITION==m);
# length_of_beginning = 11 + floor(log10(m));
# cluster_size = length(indsInM);
# if names
# text = ['Cluster ' num2str(m) ': {' char(popnames{indsInM(1)})];
# for k = 2:cluster_size
# text = [text ', ' char(popnames{indsInM(k)})];
# end;
# else
# text = ['Cluster ' num2str(m) ': {' num2str(indsInM(1))];
# for k = 2:cluster_size
# text = [text ', ' num2str(indsInM(k))];
# end;
# end
# text = [text '}'];
# while length(text)>58
# %Take one line and display it.
# new_line = takeLine(text,58);
# text = text(length(new_line)+1:end);
# disp(new_line);
# if (fid ~= -1)
# fprintf(fid,'%s \n',[new_line]);
# fprintf(fid,'\n');
# end
# if length(text)>0
# text = [blanks(length_of_beginning) text];
# else
# text = [];
# end;
# end;
# if ~isempty(text)
# disp(text);
# if (fid ~= -1)
# fprintf(fid,'%s \n',[text]);
# fprintf(fid,'\n');
# end
# end;
# end
# if npops > 1
# disp(' ');
# disp(' ');
# disp('Changes in log(marginal likelihood) if indvidual i is moved to group j:');
# if (fid ~= -1)
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', ['Changes in log(marginal likelihood) if indvidual i is moved to group j:']); fprintf(fid, '\n');
# end
# if names
# nameSizes = zeros(ninds,1);
# for i = 1:ninds
# nimi = char(popnames{i});
# nameSizes(i) = length(nimi);
# end
# maxSize = max(nameSizes);
# maxSize = max(maxSize, 5);
# erotus = maxSize - 5;
# alku = blanks(erotus);
# ekarivi = [alku ' ind' blanks(6+erotus)];
# else
# ekarivi = ' ind ';
# end
# for i = 1:cluster_count
# ekarivi = [ekarivi ownNum2Str(i) blanks(8-floor(log10(i)))];
# end
# disp(ekarivi);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [ekarivi]); fprintf(fid, '\n');
# end
# %ninds = size(data,1)/rowsFromInd;
# changesInLogml = LOGDIFF';
# for ind = 1:ninds
# %[muutokset, diffInCounts] = laskeMuutokset(ind, rowsFromInd, data, ...
# % adjprior, priorTerm);
# %changesInLogml(:,ind) = muutokset;
# muutokset = changesInLogml(:,ind);
# if names
# nimi = char(popnames{ind});
# rivi = [blanks(maxSize - length(nimi)) nimi ':'];
# else
# rivi = [blanks(4-floor(log10(ind))) ownNum2Str(ind) ':'];
# end
# for j = 1:npops
# rivi = [rivi ' ' logml2String(omaRound(muutokset(j)))];
# end
# disp(rivi);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [rivi]); fprintf(fid, '\n');
# end
# end
# disp(' '); disp(' ');
# disp('KL-divergence matrix in PHYLIP format:');
# dist_mat = zeros(npops, npops);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [' ']); %fprintf(fid, '\n');
# fprintf(fid, '%s \n', [' ']); %fprintf(fid, '\n');
# fprintf(fid, '%s \n', ['KL-divergence matrix in PHYLIP format:']); %fprintf(fid, '\n');
# end
# maxnoalle = size(COUNTS,1);
# nloci = size(COUNTS,2);
# d = zeros(maxnoalle, nloci, npops);
# prior = adjprior;
# prior(find(prior==1))=0;
# nollia = find(all(prior==0)); %Lokukset, joissa oli havaittu vain yht?alleelia.
# prior(1,nollia)=1;
# for pop1 = 1:npops
# d(:,:,pop1) = (squeeze(COUNTS(:,:,pop1))+prior) ./ repmat(sum(squeeze(COUNTS(:,:,pop1))+prior),maxnoalle,1);
# %dist1(pop1) = (squeeze(COUNTS(:,:,pop1))+adjprior) ./ repmat((SUMCOUNTS(pop1,:)+adjprior), maxnoalle, 1);
# end
# % ekarivi = blanks(7);
# % for pop = 1:npops
# % ekarivi = [ekarivi num2str(pop) blanks(7-floor(log10(pop)))];
# % end
# ekarivi = num2str(npops);
# disp(ekarivi);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [ekarivi]); %fprintf(fid, '\n');
# end
# for pop1 = 1:npops
# % rivi = [blanks(2-floor(log10(pop1))) num2str(pop1) ' '];
# for pop2 = 1:pop1-1
# dist1 = d(:,:,pop1); dist2 = d(:,:,pop2);
# div12 = sum(sum(dist1.*log2((dist1+10^-10) ./ (dist2+10^-10))))/nloci;
# div21 = sum(sum(dist2.*log2((dist2+10^-10) ./ (dist1+10^-10))))/nloci;
# div = (div12+div21)/2;
# % rivi = [rivi kldiv2str(div) ' '];
# dist_mat(pop1,pop2) = div;
# end
# % disp(rivi);
# % if (fid ~= -1)
# % fprintf(fid, '%s \n', [rivi]); fprintf(fid, '\n');
# % end
# end
# dist_mat = dist_mat + dist_mat'; % make it symmetric
# for pop1 = 1:npops
# rivi = ['Cluster_' num2str(pop1) ' '];
# for pop2 = 1:npops
# rivi = [rivi kldiv2str(dist_mat(pop1,pop2)) ' '];
# end
# disp(rivi);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [rivi]); %fprintf(fid, '\n');
# end
# end
# end
# disp(' ');
# disp(' ');
# disp('List of sizes of 10 best visited partitions and corresponding log(ml) values');
# if (fid ~= -1)
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', ['List of sizes of 10 best visited partitions and corresponding log(ml) values']); fprintf(fid, '\n');
# end
# partitionSummary = sortrows(partitionSummary,2);
# partitionSummary = partitionSummary(size(partitionSummary,1):-1:1 , :);
# partitionSummary = partitionSummary(find(partitionSummary(:,2)>-1e49),:);
# if size(partitionSummary,1)>10
# vikaPartitio = 10;
# else
# vikaPartitio = size(partitionSummary,1);
# end
# for part = 1:vikaPartitio
# line = [num2str(partitionSummary(part,1)) ' ' num2str(partitionSummary(part,2))];
# disp(line);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [line]); fprintf(fid, '\n');
# end
# end
# if ~fixedK
# disp(' ');
# disp(' ');
# disp('Probabilities for number of clusters');
# if (fid ~= -1)
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', [' ']); fprintf(fid, '\n');
# fprintf(fid, '%s \n', ['Probabilities for number of clusters']); fprintf(fid, '\n');
# end
# npopsTaulu = unique(partitionSummary(:,1));
# len = length(npopsTaulu);
# probs = zeros(len,1);
# partitionSummary(:,2) = partitionSummary(:,2)-max(partitionSummary(:,2));
# sumtn = sum(exp(partitionSummary(:,2)));
# for i=1:len
# npopstn = sum(exp(partitionSummary(find(partitionSummary(:,1)==npopsTaulu(i)),2)));
# probs(i) = npopstn / sumtn;
# end
# for i=1:len
# if probs(i)>1e-5
# line = [num2str(npopsTaulu(i)) ' ' num2str(probs(i))];
# disp(line);
# if (fid ~= -1)
# fprintf(fid, '%s \n', [line]); fprintf(fid, '\n');
# end
# end
# end
# end
# if (fid ~= -1)
# fclose(fid);
# else
# diary off
# end
# %--------------------------------------------------------------- # %---------------------------------------------------------------
@ -1488,49 +1208,6 @@ greedyMix <- function(
# num = round(num); # num = round(num);
# num2 = num/10; # num2 = num/10;
# %---------------------------------------------------------
# function digit = palautaYks(num,yks)
# % palauttaa luvun num 10^yks termin kertoimen
# % string:in?
# % yks t<>ytyy olla kokonaisluku, joka on
# % v<>hint<6E><74>n -1:n suuruinen. Pienemmill?
# % luvuilla tapahtuu jokin py<70>ristysvirhe.
# if yks>=0
# digit = rem(num, 10^(yks+1));
# digit = floor(digit/(10^yks));
# else
# digit = num*10;
# digit = floor(rem(digit,10));
# end
# digit = num2str(digit);
# function mjono = kldiv2str(div)
# mjono = ' ';
# if abs(div)<100
# %Ei tarvita e-muotoa
# mjono(6) = num2str(rem(floor(div*1000),10));
# mjono(5) = num2str(rem(floor(div*100),10));
# mjono(4) = num2str(rem(floor(div*10),10));
# mjono(3) = '.';
# mjono(2) = num2str(rem(floor(div),10));
# arvo = rem(floor(div/10),10);
# if arvo>0
# mjono(1) = num2str(arvo);
# end
# else
# suurinYks = floor(log10(div));
# mjono(6) = num2str(suurinYks);
# mjono(5) = 'e';
# mjono(4) = palautaYks(abs(div),suurinYks-1);
# mjono(3) = '.';
# mjono(2) = palautaYks(abs(div),suurinYks);
# end
# %-------------------------------------------------------------------- # %--------------------------------------------------------------------
# function dist2 = laskeOsaDist(inds2, dist, ninds) # function dist2 = laskeOsaDist(inds2, dist, ninds)