From 4a1d6d0d61894abc8b2e65e02599faf0a567d3df Mon Sep 17 00:00:00 2001 From: hamazaki1990 Date: Fri, 9 Jun 2017 11:15:48 +0900 Subject: [PATCH] solved #22 --- hamazaki1990/individual.py | 2 +- hamazaki1990/population.py | 8 ++--- hamazaki1990/visualize.py | 61 ++++++++++++++++++++------------------ hamazaki1990/visualize.r | 17 ++++++----- 4 files changed, 47 insertions(+), 41 deletions(-) diff --git a/hamazaki1990/individual.py b/hamazaki1990/individual.py index 1d42b0b..a06c778 100644 --- a/hamazaki1990/individual.py +++ b/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 diff --git a/hamazaki1990/population.py b/hamazaki1990/population.py index b4ab22e..bae86e7 100644 --- a/hamazaki1990/population.py +++ b/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() diff --git a/hamazaki1990/visualize.py b/hamazaki1990/visualize.py index 44f9946..7857cf3 100644 --- a/hamazaki1990/visualize.py +++ b/hamazaki1990/visualize.py @@ -2,37 +2,38 @@ 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 @@ -40,14 +41,16 @@ 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) diff --git a/hamazaki1990/visualize.r b/hamazaki1990/visualize.r index 9e100f6..fdf9975 100644 --- a/hamazaki1990/visualize.r +++ b/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)