using System; namespace Science.Biology.PopulationGenetics { /// /// ChromosomePair /// public class ChromosomePair { private int numberOfAllele; private int[,] configuration; // NumberOfAllele X 2 private double[,] mutation; // mutation between chromosomePairs private int[,] allele; // chromosome.NumberOfAllele X chromosome.NumberOfAllele public int NumberOfAllele { get{return numberOfAllele;} } public int[,] Configuration // chromosome configuration { get{return configuration;} } public double[,] Mutation // property of mutation { get{return mutation;} } public int[,] Allele // property of allele { get{return allele;} } private Chromosome chromosome; public ChromosomePair(Chromosome inputChromosome) { chromosome = inputChromosome; numberOfAllele = (int) Math.Floor((chromosome.NumberOfAllele+1.0)*chromosome.NumberOfAllele/2.0+0.1); FindConfiguration(); FindAllele(); FindMutation(); } private void FindConfiguration() { configuration = new int[numberOfAllele,2]; int k = 0; for(int i = 0; i < chromosome.NumberOfAllele; i++) { for(int j = i; j < chromosome.NumberOfAllele; j++) { configuration[k,0] = i; configuration[k++,1] = j; } } } private void FindAllele() { allele = new int[chromosome.NumberOfAllele,chromosome.NumberOfAllele]; int k = 0; for(int i = 0; i < chromosome.NumberOfAllele; i++) { for(int j = i; j < chromosome.NumberOfAllele; j++) { allele[i,j] = k++; allele[j,i] = allele[i,j]; } } } private void FindMutation() { mutation = new Double[numberOfAllele,numberOfAllele]; double sum; double[,] mu = new Double[numberOfAllele,numberOfAllele]; double[,] re = new Double[numberOfAllele,numberOfAllele]; Findmu(mu); if(chromosome.NumberOfGenes <= 1) TriviallyFindre(re); else Findre(re); for(int j1 = 0; j1 < numberOfAllele; j1++) { for(int j2 = 0; j2 < numberOfAllele; j2++) { sum = 0.0; for(int j3 = 0; j3 < numberOfAllele; j3++) sum += mu[j1,j3]*re[j3,j2]+re[j1,j3]*mu[j3,j2]; mutation[j1,j2] = sum/2.0; } } } private void Findmu(double[,] mu) { int i1,i2,i3,i4; for(int from = 0; from < numberOfAllele; from++) { for(int to = 0; to < numberOfAllele; to++) { i1 = configuration[from, 0]; i2 = configuration[from, 1]; i3 = configuration[to, 0]; i4 = configuration[to, 1]; if(i3 == i4) mu[from,to] = chromosome.Mutation[i1,i3]*chromosome.Mutation[i2,i4]; else mu[from,to] = chromosome.Mutation[i1,i3]*chromosome.Mutation[i2,i4] +chromosome.Mutation[i1,i4]*chromosome.Mutation[i2,i3]; } } } private void TriviallyFindre(double[,] re) { for(int i = 0; i < numberOfAllele; i++) re[i,i] = 1.0; } private void Findre(double[,] re) { int temp; double[,,] retemp = new double[numberOfAllele,numberOfAllele,chromosome.NumberOfGenes - 1]; int[] gcl = new int[chromosome.NumberOfGenes]; int[] gcr = new int[chromosome.NumberOfGenes]; for(int where = 0; where < chromosome.NumberOfGenes - 1; where++) { for(int from = 0; from < numberOfAllele; from++) { for(int to = 0; to < numberOfAllele; to++) { for(int k = 0; k < chromosome.NumberOfGenes; k++) { gcl[k] = chromosome.Configuration[configuration[from, 0],k]; gcr[k] = chromosome.Configuration[configuration[from, 1],k]; } for(int p = 0; p <= where; p++) { temp = gcl[p]; gcl[p] = gcr[p]; gcr[p] = temp; } if(to == allele[chromosome.Allele(gcl),chromosome.Allele(gcr)]) retemp[from,to,where] = chromosome.Recombination[where]; else retemp[from,to,where] = 0.0; } retemp[from,from,where] += 1.0 - chromosome.Recombination[where]; } } double[,] prod = new double[numberOfAllele,numberOfAllele]; for(int where = chromosome.NumberOfGenes - 2; where >= 1; where--) { double sum; for(int j1 = 0; j1 < numberOfAllele; j1++) { for(int j2 = 0; j2 < numberOfAllele; j2++) { sum = 0.0; for(int j3 = 0; j3 < numberOfAllele; j3++) sum += retemp[j1,j3,where-1]*retemp[j3,j2,where]; prod[j1,j2] = sum; } } for(int j1 = 0; j1 < numberOfAllele; j1++) { for(int j2 = 0; j2 < numberOfAllele; j2++) { retemp[j1,j2,where-1] = prod[j1,j2]; } } } for(int j1 = 0; j1 < numberOfAllele; j1++) { for(int j2 = 0; j2 < numberOfAllele; j2++) { re[j1,j2] = retemp[j1,j2,0]; } } } } }