solved heavywatal#22

hamazaki1990/individual.py

@@ -3,7 +3,7 @@
 
 
 class Individual:
-    mutationrate = 0.1
+    mutationrate = 0.01
 
     def __init__(self, n, f=1.0):
         self._id = n

hamazaki1990/population.py

@@ -61,7 +61,7 @@ def list_mutation(self):
                 m_list[i][k] += 1
         return m_list
 
-    def calculate_mutantfreq_per_site(self):
+    def mutantfreq_per_site(self):
         genotypes = [x.get_genotype() for x in self._inds]
         m_sites = []
         for x in genotypes:
@@ -73,8 +73,8 @@ def calculate_mutantfreq_per_site(self):
                 k = m_sites.index(genotypes[i][j])
                 m_count[k] += 1
         mutantfreq = [x/len(self._inds) for x in m_count]
-        mutantfreq_per_site = ({m_sites[i]: mutantfreq[i]
-                               for i in range(len(m_sites))})
+        mutantfreq_per_site = ([[m_sites[i], mutantfreq[i]]
+                               for i in range(len(m_sites))])
         return mutantfreq_per_site
 
     def is_not_fixed(self):
@@ -96,7 +96,7 @@ def main():
         print(p1_1.get_ids())
         print(p1_1.get_genotypes())
         print(p1_1.list_mutation())
-        print(p1_1.calculate_mutantfreq_per_site())
+        print(p1_1.mutantfreq_per_site())
 
     print(p1_1.is_not_fixed())
     p1_1.list_mutation()

hamazaki1990/visualize.py

@@ -2,52 +2,55 @@
 from population import Population
 
 
+def get_step(t, mutations):
+    t_mutations = []
+    for i in mutations:
+        t_mutation = [t]
+        t_mutation.extend(i)
+        t_mutations.append(t_mutation)
+    return t_mutations
+
+
 def change_allelewf(population, generation):
-    change_allele = [population.calculate_mutantfreq_per_site()]
+    t = 0
+    change_allele = population.mutantfreq_per_site()
+    derived_allelefreq = get_step(t, change_allele)
     for x in range(generation):
+        t += 1
         population.next_genwf()
-        change_allele.append(population.calculate_mutantfreq_per_site())
-    all_mutations = []
-    for x in change_allele:
-        all_mutations.extend(x.keys())
-    all_mutationsites = sorted(list(set(all_mutations)))
-    derived_allelefreq = []
-    for x in range(generation+1):
-        t = [x]
-        t.extend([change_allele[x].get(y, 0.0) for y in all_mutationsites])
-        derived_allelefreq.append(t)
+        change_allele = population.mutantfreq_per_site()
+        next_gen = get_step(t, change_allele)
+        derived_allelefreq.extend(next_gen)
     return derived_allelefreq
 
 
 def change_allelemo(population, generation):
-    change_allele = [population.calculate_mutantfreq_per_site()]
+    t = 0
+    change_allele = population.mutantfreq_per_site()
+    derived_allelefreq = get_step(t, change_allele)
     for x in range(generation):
+        t += 1
         population.next_genmo()
-        change_allele.append(population.calculate_mutantfreq_per_site())
-    all_mutations = []
-    for x in change_allele:
-        all_mutations.extend(x.keys())
-    all_mutationsites = sorted(list(set(all_mutations)))
-    derived_allelefreq = []
-    for x in range(generation+1):
-        t = [x]
-        t.extend([change_allele[x].get(y, 0.0) for y in all_mutationsites])
-        derived_allelefreq.append(t)
+        change_allele = population.mutantfreq_per_site()
+        next_gen = get_step(t, change_allele)
+        derived_allelefreq.extend(next_gen)
     return derived_allelefreq
 
 
 def output_allelechange(filename, function, population, generation):
     with open(filename, "w", encoding="utf-8") as outfile:
         writer = csv.writer(outfile)
         freq = function(population, generation)
-        num = ["allele" + str(i) + "frequency" for i in range(len(freq[0])-1)]
-        col_name = ["generation"]
-        col_name.extend(num)
+        col_name = ["step", "locus", "frequency"]
         writer.writerow(col_name)
         for x in range(len(freq)):
             writer.writerow(freq[x])
 
 
+example = [[0.1, 0.2], [0.2, 0.2], [0.8, 0.1]]
+print(example[1])
+print(get_step(2, example))
+
 p1_1 = Population(10)
 fixprocess1 = change_allelewf(p1_1, 10)
 print(fixprocess1)
@@ -56,8 +59,8 @@ def output_allelechange(filename, function, population, generation):
 fixprocess2 = change_allelemo(p1_2, 20)
 print(fixprocess2)
 
-p1 = Population(100)
-output_allelechange("allelefreqwf.csv", change_allelewf, p1, 100)
+p1 = Population(1000)
+output_allelechange("allelefreqwf.csv", change_allelewf, p1, 1000)
 
-p2 = Population(100)
-output_allelechange("allelefreqmo.csv", change_allelemo, p2, 100)
+p2 = Population(1000)
+output_allelechange("allelefreqmo.csv", change_allelemo, p2, 1000)

hamazaki1990/visualize.r

@@ -1,13 +1,16 @@
 library("tidyverse")
 setwd("~/github/practice-py/hamazaki1990/")
 freq_wf <- read_csv("allelefreqwf.csv") %>%
-    dplyr::mutate(model = "wf")
+  dplyr::arrange(locus, step)
 freq_wf
+plotwf <- ggplot(freq_wf, aes(x=step, y=frequency, group=locus))
+plotwf <-plotwf +geom_line()
+plotwf <- plotwf + ylim(0, 1)
+print(plotwf)
 freq_mo <- read_csv("allelefreqmo.csv") %>%
-    dplyr::mutate(model = "moran")
+  dplyr::arrange(locus, step)
 freq_mo
-freq <- dplyr::bind_rows(freq_wf, freq_mo)
-print(freq)
-freqplot <- ggplot(freq, aes(x=generation, y=derived_allele_frequency, group=model, color=model))
-freqplot <-freqplot + geom_line()
-print(freqplot)
+plotmo <- ggplot(freq_mo, aes(x=step, y=frequency, group=locus))
+plotmo <-plotmo +geom_line()
+plotmo <- plotmo + ylim(0, 1)
+print(plotmo)