Skip to content

Instantly share code, notes, and snippets.

@jeanpat

jeanpat/Cellular.py

Last active July 15, 2018 16:50
Show Gist options
  • Save jeanpat/472053f4c12bd06e214f5944b49d5f65 to your computer and use it in GitHub Desktop.
Save jeanpat/472053f4c12bd06e214f5944b49d5f65 to your computer and use it in GitHub Desktop.
Download ZIP
Modules used in telomere length modelisation
Raw
Cellular.py
# -*- coding: utf-8 -*-
"""
Created on Wed Mar 30 13:04:37 2016
@author: jeanpat
"""
import pandas as pn
from itertools import product
import ClassChromosome as K
import Genome as G
def init_genome(telomeres_distribution =
{"1":(15000,18000,9000,12000),
"2":(11000, 20000, 11050,19950)}):
'''telomeres_distribution is a dictionary:
{"1":(15000,18000,9000,12000), "2":(11000, 20000, 11050,19950)}
'''
if type(telomeres_distribution) != type(dict()):
raise NameError('Telomere length must be a dictionary "chrom1":(p1,q1,p2,q2)')
chrom_set = []
for tels in telomeres_distribution.keys():
t0 = telomeres_distribution[tels][0]#first homolog at pter
t1 = telomeres_distribution[tels][1]#first homolog at qter
t2 = telomeres_distribution[tels][2]#other homolog at pter
t3 = telomeres_distribution[tels][3]#other homolog at qter
w1 = K.Watson(ttaggg = t0, ccctaa = t1)
c1 = K.Crick(ttaggg = t1, ccctaa = t0)
w2 = K.Watson(ttaggg= t2, ccctaa = t3)
c2 = K.Crick(ttaggg= t3, ccctaa = t2)
#chrom_set.ap
chrom_set.append(K.Chromosome(watson=w1,crick=c1, name =str(tels),homolog='mat'))
chrom_set.append(K.Chromosome(watson=w2,crick=c2, name =str(tels),homolog='pat'))
return G.Genome(*chrom_set)
class Cell(object):
# def init_genome(self, telomeres_distribution = {"1":(15000,18000,9000,12000),
# "2":(11000, 20000, 11050,19950)}):
# '''telomeres_distribution is a dictionary:
# {"1":(15000,18000,9000,12000), "2":(11000, 20000, 11050,19950)}
#
# '''
# print "init_genome"
# if type(telomeres_distribution)<>type(dict()):
# raise NameError('Telomere length must be a dictionary "chrom1":(p1,q1,p2,q2)')
# chrom_set = []
# for tels in telomeres_distribution.keys():
# t0 = telomeres_distribution[tels][0]#first homolog at pter
# t1 = telomeres_distribution[tels][1]#first homolog at qter
# t2 = telomeres_distribution[tels][2]#other homolog at pter
# t3 = telomeres_distribution[tels][3]#other homolog at qter
#
# w1 = K.Watson(ttaggg = t0, ccctaa = t1)
# c1 = K.Crick(ttaggg = t1, ccctaa = t0)
#
# w2 = K.Watson(ttaggg= t2, ccctaa = t3)
# c2 = K.Crick(ttaggg= t3, ccctaa = t2)
# #chrom_set.ap
# chrom_set.append(K.Chromosome(watson=w1,crick=c1, name =str(tels),homolog='mat'))
# chrom_set.append(K.Chromosome(watson=w2,crick=c2, name =str(tels),homolog='pat'))
# return G.Genome(*chrom_set)
#
def __init__(self,division = 0,
genome = init_genome(),
short_telomere_threshold = 200,
telomerase = {'ON':False,
'Docking':(1000,0.95),
'Elongation':(200,0.95)}
):
self.division = division
#
#telomerase parametres
self.telomerase_par = telomerase
self.telomerase_ON = self.telomerase_par['ON']
#telomere lenght theshold (bp)for telomerase docking
#simulation of the D-loop opening:
self.docking_threshold = self.telomerase_par['Docking'][0]
#docking probability: bellow the docking threshold,
#the probability for the telomerase to dock is:
self.docking_proba = self.telomerase_par['Docking'][1]
# When telomerase is docked the TTAGGG motif is elongated
# by n bases (ex:200) with a probability of p (p=0.95):
self.tel_elongation = self.telomerase_par['Elongation'][0]
self.elongation_prob = self.telomerase_par['Elongation'][1]
#### GENOME initialization ####
self.genome = genome
self.short_telomere_threshold = short_telomere_threshold
self.genome_dataframe()
self.shortest_telo = self.scan_shortest_telomere()
if self.shortest_telo <= self.short_telomere_threshold:
## The cell becomes senescent, blocked in G0
self.G0 = True
else:
self.G0 = False
#def scan_shortest_telomere(self):
# df = self.genome_dataframe()
# return df.T.min()[0]
def scan_shortest_telomere(self):
df = self.genome_dataframe()
#telomeres_list = df.columns-['division']
telomeres_list = df.columns[1:]
shortest_telo = df.ix[:, telomeres_list ].min().min()
#print shortest_telo
return shortest_telo
def mitose(self):
# check if G_0 can be left
#print type(self.genome)
# trigger genome replication
self.shortest_telo = self.scan_shortest_telomere()
#print type(shortest_telo[0]),shortest_telo[0]
if self.shortest_telo <= self.short_telomere_threshold:
## The cell becomes senescent, blocked in G0
self.G0 = True
#return None
else:
self.genome.set_G1_to_G2()
#trigger telophase
son_genome = self.genome.G2_to_M()
#print self.division
self.division = self.division +1
self.genome_dataframe()
daughter_cell = Cell(division = self.division,
genome = son_genome,
telomerase = self.telomerase_par,
short_telomere_threshold = self.short_telomere_threshold)
# Check post replication if the shortest telomere is bellow threshold
#Check mother cell
self.shortest_telo = self.scan_shortest_telomere()
if self.shortest_telo <= self.short_telomere_threshold:
## The cell becomes senescent, blocked in G0
self.G0 = True
#Daughter cell: the telomere length is checked on initialisation. If the shortest telo
#of the daughter cell is bellow threshold, the daughter cell is blocked in G0
return daughter_cell
def genome_dataframe(self):
'''Return all the telomere in a pandas data_frame, one column per telomere
the homologs must be distinguished by a tag: eg:maternal,paternal
'''
chroms = self.genome.get_chromosomes_list()
homogs = self.genome.get_homologous_tags()
telo = {}
#telo['division'] = self.division
div_df = pn.DataFrame({'division' : self.division}, index = [0])
for telomeric in product(chroms,['pter','qter'], '_',homogs):
telo[''.join(telomeric)]= [self.genome.get_telo_length(pair=telomeric[0],
hog=telomeric[3],
tel=telomeric[1])]
self.telo_df = pn.DataFrame(telo)#pandas dataframe
self.telo_df = pn.concat([div_df,self.telo_df], axis = 1)
return self.telo_df
Raw
ClassChromosome.py
# -*- coding: utf-8 -*-
"""
Created on Thu Sep 6 11:29:52 2012
@author: Jean-Patrick Pommier
"""
from WatsonCrick import Watson, Crick
import random
class Chromosome(object):
def __init__(self,watson = Watson(ttaggg = 10001, ccctaa = 5000), crick = Crick(ttaggg = 5000, ccctaa = 10000), name ='1', homolog ='a'):
#print 'new chromosome'
watson.in_use = True
crick.in_use = True
self.name = name
self.homolog = homolog
self.phase = 'G1'
self.watson = watson
self.crick = crick
#Double strand DNA
self.chromatide = []
self.chromatide.append(self.watson)
self.chromatide.append(self.crick)
self.chromosome = []
self.chromosome.append(self.chromatide)#chromosome with one chromatid
def set_G1_to_G2(self):
self.phase = 'G2'
self.sister_chromatid1 = []
self.sister_chromatid2 = []
# temporary variables
self.ctd = self.chromosome.pop()
self.strand1 = self.ctd.pop()
self.strand2 = self.ctd.pop()
self.sister_chromatid1.append(self.strand1)
self.sister_chromatid1.append(self.strand1.replicate())
self.sister_chromatid2.append(self.strand2)
self.sister_chromatid2.append(self.strand2.replicate())
self.chromosome.append(self.sister_chromatid1)
self.chromosome.append(self.sister_chromatid2)
def segregate(self) :
'''gives two sister chromatids from metaphasic chromosome
'''
def getWatsonCrick_Index_in(chromatid):
if type(chromatid[0]) == Watson:
return 0, 1
else:
return 1, 0
if self.phase == 'G2':
index = random.randint(0,1)
#print index
# get one of the two chromatids of the G2 chromosome
self.sister_ctd1 = self.chromosome.pop(index)
w, c = getWatsonCrick_Index_in(self.sister_ctd1)
#Back to G1 state in the cell cycle after chromosome segregation
self.phase = 'G1'
son = Chromosome(self.sister_ctd1.pop(w),self.sister_ctd1.pop(), name = self.name, homolog=self.homolog)
return son
else:
raise NameError ('chromosome segregation only on metaphasic chrom!')
#==============================================================================
#
#==============================================================================
def get_Crick(self):
'''Returns a single strand of type Crick
'''
if type(self.chromosome[0][0]) == Crick:
return self.chromosome[0][0]
else:
return self.chromosome[0][1]
def get_Watson(self):
'''Returns a single strand of type Watson
'''
if type(self.chromosome[0][0]) == Watson:
return self.chromosome[0][0]
else:
return self.chromosome[0][1]
def get_pter_G1(self, motif='ccctaa'):
'''telomere length at pter = nber of ccctaa motifs on Crick single strand
telomere length at pter = nber of ttaggg motifs on Watson single strand
'''
#get index of crick strand in chromatid[]
if motif == 'ccctaa':
if self.phase == 'G1':
if type(self.chromosome[0][0]) == Crick:
return self.chromosome[0][0].getCCCTAA()
else:
return self.chromosome[0][1].getCCCTAA()
else:
raise NameError('only on G1 chromosome')
elif motif =='ttaggg':
if self.phase == 'G1':
if type(self.chromosome[0][0]) == Watson:
return self.chromosome[0][0].getTTAGGG()
else:
return self.chromosome[0][1].getTTAGGG()
else:
raise NameError('only on G1 chromosome')
def get_qter_G1(self,motif='ccctaa'):
'''telomere length at qter = nber of ccctaa motifs on Watson single strand
telomere length at qter = nber of ttaggg motifs on Crick single strand
'''
#get index of crick strand in chromatid[]
if motif == 'ccctaa':
if self.phase == 'G1':
if type(self.chromosome[0][0]) == Watson:
return self.chromosome[0][0].getCCCTAA()
else:
return self.chromosome[0][1].getCCCTAA()
else:
raise NameError('only on G1 chromosome')
elif motif =='ttaggg':#ERROR
if self.phase == 'G1':
if type(self.chromosome[0][0]) == Crick:
return self.chromosome[0][0].getTTAGGG()
else:
return self.chromosome[0][1].getTTAGGG()
else:
raise NameError('only on G1 chromosome')
def get_pter_G2(self,motif='ccctaa'):
'''telomere length at pter = nber of ccctaa motifs on Crick single strand
telomere length at pter = nber of ttaggg motifs on Watson single strand
'''
def getCrickindex_chromatid1():
'''get the index of the Crick strand in chromatid1
'''
if type(self.sister_chromatid1[0]) == Crick:
return 0
else:
return 1
def getCrickindex_chromatid2():
'''get the index of the Crick strand in chromatid2
'''
if type(self.sister_chromatid2[0]) == Crick:
return 0
else:
return 1
def getWatsonindex_chromatid1():
if type(self.sister_chromatid1[0]) == Watson:
return 0
else:
return 1
def getWatsonindex_chromatid2():
if type(self.sister_chromatid2[0]) == Watson:
return 0
else:
return 1
if self.phase == 'G2':
if motif == 'ccctaa':
c1 = getCrickindex_chromatid1()
c2 = getCrickindex_chromatid2()
return self.sister_chromatid1[c1].getCCCTAA(),self.sister_chromatid2[c2].getCCCTAA()
elif motif == 'ttaggg':
w1 = getWatsonindex_chromatid1()
w2 = getWatsonindex_chromatid2()
ttaggg1 = self.sister_chromatid1[w1].getTTAGGG()
ttaggg2 = self.sister_chromatid2[w2].getTTAGGG()
return ttaggg1, ttaggg2
else:
raise NameError('only on G2 chromosome')
def get_qter_G2(self,motif='ccctaa'):
'''telomere length at qter = nber of ccctaa motifs on Watson single strand
telomere length at qter = nber of ttaggg motifs on Crick single strand
on each chromatid
'''
#get index of Watson strand in chromatid[]
def getWatsonindex_chromatid1():
'''get the index of the Watson strand in chromatid1
'''
if type(self.sister_chromatid1[0]) == Watson:
return 0
else:
return 1
def getWatsonindex_chromatid2():
'''get the index of the Watson strand in chromatid2
'''
if type(self.sister_chromatid2[0]) == Watson:
return 0
else:
return 1
def getCrickindex_chromatid1():
'''get the index of the Crick strand in chromatid1
'''
if type(self.sister_chromatid1[0]) == Crick:
return 0
else:
return 1
def getCrickindex_chromatid2():
'''get the index of the Crick strand in chromatid2
'''
if type(self.sister_chromatid2[0]) == Crick:
return 0
else:
return 1
if self.phase == 'G2':
if motif == 'ccctaa':
w1 = getWatsonindex_chromatid1()
w2 = getWatsonindex_chromatid2()
return self.sister_chromatid1[w1].getCCCTAA(),self.sister_chromatid2[w2].getCCCTAA()
elif motif == 'ttaggg':
c1 = getCrickindex_chromatid1()
c2 = getCrickindex_chromatid2()
ttaggg1 = self.sister_chromatid1[c1].getTTAGGG()
ttaggg2 = self.sister_chromatid2[c2].getTTAGGG()
return ttaggg1, ttaggg2
else:
raise NameError('only on G2 chromosome')
def telomere_length(self,end = 'pter', motif = 'ccctaa'):
'''returns the telomere length (TTAGGG with motif='ttaggg' or CCCTAA with motif='ccctaa')
of a chromosome, given the location:end='pter' or end='qter'.
For a chromosome in G1 cell cycle, the length is a single number.
For a chromosome in G2, the length is a tuple, one value for each chromatid.
'''
if self.phase == 'G1':
if end == 'pter':
if motif =='ccctaa':
return self.get_pter_G1(motif)
elif motif =='ttaggg':
return self.get_pter_G1(motif)
elif end =='qter':
if motif =='ccctaa':
return self.get_qter_G1(motif)
elif motif =='ttaggg':
return self.get_qter_G1(motif)
elif self.phase =='G2':
if end == 'pter':
if motif =='ccctaa':
return self.get_pter_G2(motif)
elif motif =='ttaggg':
return self.get_pter_G2(motif)
elif end =='qter':
if motif =='ccctaa':
return self.get_qter_G2(motif)
elif motif =='ttaggg':
return self.get_qter_G2(motif)
#==============================================================================
#
# S
#
#
#==============================================================================
def test_get_pter():
k = Chromosome()
print( k)
print( 'G1+ccctaa:k.get_pter_G1() ',k.get_pter_G1())
print( "G1+ccctaa:k.get_pter_G1(motif='ccctaa') ",k.get_pter_G1(motif='ccctaa'))
print( "G1+ttaggg:k.get_pter_G1(motif='ttaggg') ",k.get_pter_G1(motif='ttaggg'))
print()
print( '----------------------')
k.set_G1_to_G2()
print( k)
print( "G2+ccctaa:k.get_pter_G1(motif='ccctaa') ",k.get_pter_G2(motif='ccctaa'))
print( "G2+ttaggg:k.get_pter_G1(motif='ttaggg') ",k.get_pter_G2(motif='ttaggg'))
print()
print()
def test_one_replication():
k = Chromosome()
#print k.chromosome
print( 'telomere p/q length G1 (ccctaa):',k.get_pter_G1(motif='ccctaa'),' ',k.get_qter_G1(motif='ccctaa'))
print( 'telomere p/q length G1 (ttaggg):',k.get_pter_G1(motif='ttaggg'),' ',k.get_qter_G1(motif='ttaggg'))
k.set_G1_to_G2()
motif='ccctaa'
print( 'ccctaa G2:','pter',k.get_pter_G2(motif),'qter',k.get_qter_G2(motif))
motif='ttaggg'
print( 'ttaggg G2:','pter',k.get_pter_G2(motif),'qter',k.get_qter_G2(motif))
print( k.chromosome)
neo = k.segregate()
print( 'neo',neo.get_pter_G1(),' ',neo.get_qter_G1())
print( 'c post seg:',k.get_pter_G1(),' ',k.get_qter_G1())
def test_replication(n):
k = Chromosome()
for i in range(n):
k.set_G1_to_G2()
k1 = k.segregate()
print( i,' ',k.get_pter_G1(),',',k.get_qter_G1(),' ',k1.get_pter_G1())
def test_rep(n):
import numpy as np
klist0 = []
k = Chromosome()
print( 'chrom name', k.name)
klist0.append(k)
generation = {}
motif = 'ttaggg'
generation[0] = (k.get_pter_G1(motif), k.get_qter_G1(motif))
for i in range(n):
klist1 = []
for c in klist0:
c.set_G1_to_G2()
for c in klist0:
klist1.append(c.segregate())
klist0 = klist0 + klist1
pq = []
for c in klist0:
pq.append((c.get_pter_G1(motif), c.get_qter_G1(motif)))
generation[i+1] = pq
# print generation.keys()
# print ((0,generation[0][0],generation[0][1]))
table = np.asarray((0,generation[0][0],generation[0][1]))
# print table.dtype
# print 'gen 0:',generation[0]
for g in generation.keys():
# print g,type(generation[g])
for pair in generation[g]:
if g == 0:
break
else:
# print 'g:',g,'pair:',pair,' ',pair[0],pair[1]
tmp = np.asarray((g,pair[0],pair[1]))
table = np.vstack([table,tmp])
return table
#==============================================================================
# Tests
#==============================================================================
if __name__ == "__main__":
test_get_pter()
test_one_replication()
telomeres=test_rep(9)
pter=telomeres[:,1]
print( pter)
qter=telomeres[:,2]
gen=telomeres[:,0]
import matplotlib.pyplot as plt
#import matplotlib.mlab as mlab
from scipy import stats
#import scipy as sc
fig = plt.figure(1)
plt.subplot(312, frameon = True)
plt.scatter(pter,qter,c=gen,s=50)
slope, intercept, r_value, p_value, slope_std_error = stats.linregress(gen, pter)
print( 'slope, intercept:',slope, intercept)
print( telomeres[gen ==4,1])
#plt.scatter(gen,pter)
plt.subplot(311, frameon = True)
# ax = fig.add_subplot(311)
# mu = sc.mean
# y = mlab.normpdf(bincenters, mu, sigma)
plt.hist(telomeres[gen==9,1], bins=100)
plt.subplot(313, frameon = True)
plt.scatter(gen, pter)
plt.show()
Raw
Genome.py
# -*- coding: utf-8 -*-
"""
Created on Thu Sep 13 18:25:58 2012
@author: JeanPat
"""
#from WatsonCrick import *
import ClassChromosome as K
import numpy as np
#import pdb
class Genome(object):
''' a genome is a collection of chromosomes
'''
def __init__(self, *args):
def tuplesChrom_to_dict(tu):
''' convert a tuple (('1','a'), ('1','A'), ('2','b'), ('2','B'))
into a dict {'1':('a','A'),'2':('b','B')}
'''
g={}
for c,krom in enumerate(tu):
if krom.name in g:
g[krom.name]=g[krom.name],krom
else:
g[krom.name]=krom
return g
#print type(args)
if len(args) % 2 ==0:
self.ploidy = (len(args) / 2)
else:
self.ploidy = len(args) / 2 , 1
#print('WARNING: Non diploïd genome, delta p/q will crash')
#========================================================================
# Let's populate the genome with chromosomes
#========================================================================
for count, chrom in enumerate(args):
#pdb.set_trace()
if type(chrom) != type(K.Chromosome()):
raise NameError('Genome only loves Chromosomes')
#==============================================================================
# self.genome is dictionary where :
# key= chromosome name (1, 2, 3...)
# values: pair of homologous chromosomes (a or b for ex)
#==============================================================================
self.genome = tuplesChrom_to_dict(args)
#==============================================================================
#
#Let's define what we can get from a genome object:
# * telomere length :
# -of a given chromosome/homolog
# -mean telomere length for the whole genome
# -std dev
# * difference of length:
# -between two homologs
# -distribution of telomere length
# * telomere length ratio:
# -pter / qter for a given chromosome/homolog
# -ratio distribution
# * ordered list of telomeres : ex 1pa, 2pB, 2qA, 1pA, 1qB ...
#==============================================================================
def get_chrom_phase(self, pair='1', hog='a'):
'''get if a chrom '1' homolog 'a' is G1 or G2
'''
homolog = self.get_homologous_chromosome(pair,hog)
return homolog.phase
def get_chromosomes_list(self):
return self.genome.keys()
def get_chrom_pair(self,pair='1'):
return self.genome[pair]#should be a tuple
def get_homologous_tags(self, pair='1'):
'''given a pair of chromosomes, returns the values
of Chromosome.homolog ex ['a','A']
'''
pair = self.get_chrom_pair(pair)
tags=[]
for k in pair:
tags.append(k.homolog)#k should be a Chromosome instance
#print 'tags ', tags
return tags
def get_homologous_chromosome(self,pair,hog):
'''Return one of the two chromosomes: ex
Need chromosome 1, homolog a? ask get_homologous_chromosome('1','a')
'''
pair = self.get_chrom_pair(pair)
chrom1 = pair[0]
chrom2 = pair[1]
if chrom1.homolog == hog:
return chrom1
elif chrom2.homolog == hog:
return chrom2
def get_telo_length(self, pair='1',hog='a', tel='pter',cycle='G1',motif='ccctaa'):
#==============================================================================
# Should use .homolog.phase to know if chromosome is in G1 or G2
#==============================================================================
homolog = self.get_homologous_chromosome(pair, hog)
if (tel == 'pter') and (cycle == 'G1'):
return homolog.get_pter_G1(motif)
elif tel == 'qter' and cycle == 'G1':
return homolog.get_qter_G1(motif)
elif tel == 'pter' and cycle == 'G2':
return homolog.get_pter_G2(motif)
elif tel == 'qter' and cycle == 'G2':
return homolog.get_qter_G2(motif)
else:
raise NameError('Can only get length at pter or qter in G1 or G2 cell cycle')
def get_telomere_length(self,pair='1',hog='a', tel='pter', motif='ccctaa'):
'''need not to know if G1 or G2
'''
homolog = self.get_homologous_chromosome(pair, hog)
return homolog.telomere_length(end = tel, motif = motif)
def get_mean_std_telo_lenght(self, cycle='G1',motif='ccctaa'):
'''compute the mean and standard deviation telomere length in one genome
'''
telomeres = []#an array containing all the telomeres
#print 'genome:',self.genome
for c,kpair in self.genome.items():
homolog1=kpair[0]
homolog2=kpair[1]
#==============================================================================
# Uses ccctaa motifs and chromosomes when G1 cycle is 'on'
#==============================================================================
p1=self.get_telomere_length(pair=homolog1.name, hog=homolog1.homolog,tel='pter')
q1=self.get_telomere_length(pair=homolog1.name, hog=homolog1.homolog,tel='qter')
p2=self.get_telomere_length(pair=homolog2.name, hog=homolog2.homolog,tel='pter')
q2=self.get_telomere_length(pair=homolog2.name, hog=homolog2.homolog,tel='qter')
telomeres.append(p1)
telomeres.append(q1)
telomeres.append(p2)
telomeres.append(q2)
tel = np.asarray(telomeres)
mean_telo = tel.mean()
std_telo = tel.std()
# print telomeres
# print (tel-mean_telo)/std_telo# centrée réduites
return mean_telo,std_telo
def get_pq_ratio(self,pair='1',hog='a',motif='ccctaa'):
'''Compute ratio (length pter)/(length qter) for
a given homologous chromsome
'''
p = self.get_telomere_length(pair,hog,motif,tel='pter')
q = self.get__telomere_length(pair,hog,motif,tel='qter')
if q>0:
return 1.0*p/q
else:
return None
def delta_telo(self,chromPair='1',tel='pter', cycle='G1',motif='ccctaa'):
''' calculate the difference of length between two homologous telomeres
'''
hmlog_tags = self.get_homologous_tags(chromPair)
first_hmlog = hmlog_tags[0]
second_hmlog = hmlog_tags[1]
if cycle == 'G1':
telo1 = self.get_telomere_length(chromPair, first_hmlog, tel,motif)
# print 'telo1',telo1
telo2 = self.get_telomere_length(chromPair, second_hmlog, tel,motif)
# print 'telo2',telo2
# print telo1-telo2
return telo1-telo2
elif cycle == 'G2':
pass
else:
raise NameError('only G1 or G2 cell cycle state')
def distribution_delta_telo(self, mode='raw'):
''' Should return all the difference between all the telomeres two by two
'''
pass
#==================================================
# Let's define what we can do with a genome
# * cell cycle from G1 to G2
# * segregating :1 Genome -----> 2 Genomes
#==================================================
def set_G1_to_G2(self):
for n,pair_krom in self.genome.items():
pair_krom[0].set_G1_to_G2()
pair_krom[1].set_G1_to_G2()
def G2_to_M(self):
#
chromosomes_list=[]
for n, pair_krom in self.genome.items():
rosalind = pair_krom[0].segregate()# one chromatid of one homologous chromosome
franklin = pair_krom[1].segregate()# one chromatid of the other homologous chromosome
chromosomes_list.append(rosalind)
chromosomes_list.append(franklin)
return Genome(*chromosomes_list)
Display the source blob
Display the rendered blob
Raw
Telomere Length Distribution in Synchronous Dividing Cells.ipynb
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Raw
unittest-watsoncrick.py
# -*- coding: utf-8 -*-
"""
Created on Fri Oct 5 11:05:01 2012
@author: Jean-Pat
"""
import unittest as ut
import WatsonCrick as ss# single strand DNA
class Test_Wa(ut.TestCase):
""" class for unit tests for Watson class """
def test_initW(self):
""" Watson constructor test"""
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
self.assertEqual(w.getTTAGGG(), 10400)
self.assertEqual(w.getCCCTAA(), 5000)
def test_initCr(self):
""" Crick constructor test"""
c = ss.Crick(ttaggg = 5000, ccctaa = 10400)
self.assertEqual(c.getTTAGGG(), 5000)
self.assertEqual(c.getCCCTAA(), 10400)
def test_delete5prime(self):
""" test 5' ccctaa deletion with No exonuclease"""
ss.deletion = 1
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
w.delete5prime()
self.assertLess(w.getCCCTAA(), 5000)
c = ss.Crick(ttaggg = 5000, ccctaa = 10400)
c.delete5prime()
self.assertLess(c.getCCCTAA(), 10400)
def test_replicate1(self):
""" Check that Watson replication gives Crick and reciprocally"""
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
self.assertIsInstance(w.replicate(),ss.Crick)
c = ss.Crick(ttaggg = 5000, ccctaa = 10400)
self.assertIsInstance(c.replicate(),ss.Watson)
def test_replicate2(self):
""" Check Watson.replicate(),or Crick, reduces Watson (Crick) CCCTAA"""
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
w.replicate()
self.assertLess(w.getCCCTAA(),5000)
c = ss.Crick(ttaggg = 5000, ccctaa = 10400)
c.replicate()
self.assertLess(c.getCCCTAA(), 10400)
def test_replicate3(self):
""" Take a Watson strand
replicate it
Check that the new Crick CCCTAA< Watson TTAGGG"""
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
c = w.replicate()
self.assertLess(c.getCCCTAA(),w.getTTAGGG())
def test_replicate4(self):
""" Take a Crick strand
replicate it
Check that the new Watson CCCTAA< Crick TTAGGG"""
c = ss.Crick(ttaggg = 10400, ccctaa = 5000)
w = c.replicate()
self.assertLess(w.getCCCTAA(),c.getTTAGGG())
def test_lengthenTTAGGG(self):
""" Check that Watson (Crick) TTAGGG can be lengthened (telomerase)"""
w = ss.Watson(ttaggg = 10400, ccctaa = 5000)
w.lengthen_ttaggg(200)
self.assertGreater(w.getTTAGGG(),10400)
c = ss.Crick(ttaggg = 5000, ccctaa = 10000)
c.lengthen_ttaggg(200)
self.assertGreater(c.getTTAGGG(),5000)
if __name__ == '__main__':
ut.main()
Raw
unittest_Cell.py
# -*- coding: utf-8 -*-
"""
Created on Wed Mar 30 16:52:42 2016
@author: jeanpat
"""
import unittest as ut
import Cellular as C
class test_cell(ut.TestCase):
def test_initDivision(self):
cell = C.Cell()
self.assertEqual(cell.division, 0)
def test_initDivision1(self):
cell = C.Cell()
cell.mitose()
self.assertEqual(cell.division, 1)
def test_length1(self):
tel=(15000,18000,9000,201)
distrib = {"1":tel}
tgen = C.init_genome(telomeres_distribution = distrib)
tcell = C.Cell(genome = tgen, short_telomere_threshold= 200)
df = tcell.genome_dataframe()
self.assertEqual(df['1pter_mat'][0], 15000)
def test_cell_cycle1(self):
tel=(15000,18000,9000,201)
distrib = {"1":tel}
tgen = C.init_genome(telomeres_distribution = distrib)
tcell = C.Cell(genome = tgen, short_telomere_threshold= 200)
self.assertEqual(tcell.G0, False)
def test_cell_cycle2(self):
tel=(15000,18000,9000,201)
distrib = {"1":tel}
tgen = C.init_genome(telomeres_distribution = distrib)
mother_cell = C.Cell(genome = tgen, short_telomere_threshold= 200)
daughter_cell = mother_cell.mitose()
self.assertEqual(mother_cell.G0, True)
self.assertEqual(daughter_cell.G0, True)
if __name__ == '__main__':
ut.main()
Raw
unittest_Chromosome.py
# -*- coding: utf-8 -*-
"""
Created on Fri Oct 5 15:01:23 2012
@author: Jean-Pat
"""
import unittest as ut
import ClassChromosome as K
class Test_Chromosome(ut.TestCase):
""" class for unit tests for Chromosome class """
def test_initK(self):
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k=K.Chromosome(watson=w,crick=c)
self.assertIsInstance(k,K.Chromosome)
def test_Chromatid1(self):
""" check if a chromatid (Chromosome) has two
Watson an Crick obj
"""
k=K.Chromosome()
self.assertIsInstance(k.get_Crick(),K.Crick)
self.assertIsInstance(k.get_Watson(),K.Watson)
def test_1_telomere_length(self):
""" check telomere length on instanciation
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k=K.Chromosome(w, c)
length_pter = k.get_pter_G1()
length_qter = k.get_qter_G1()
self.assertEqual(length_pter,10000)
self.assertEqual(length_qter,5000)
def test_2a_telomere_length_G2(self):
""" check telomere length after G1 to G2 is tupple
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k=K.Chromosome(w, c)
k.set_G1_to_G2()
cell_cycle ='G2'
motif='ccctaa'
self.assertEqual(k.phase,cell_cycle)
self.assertIsInstance(k.get_pter_G2(motif),tuple)
self.assertIsInstance(k.get_qter_G2(motif),tuple)
def test_2b_telomere_length_G2(self):
""" check telomere length after G1 to G2 is tupple
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k=K.Chromosome(w, c)
k.set_G1_to_G2()
cell_cycle ='G2'
motif='ttaggg'
self.assertEqual(k.phase,cell_cycle)#self.assertEquals(k.phase,cell_cycle)
self.assertIsInstance(k.get_pter_G2(motif),tuple)
self.assertIsInstance(k.get_qter_G2(motif),tuple)
def test_3_telomere_G1_pter_GStrand(self):
""" Take a chromosome :
*G1 to G2
*segregate (mitosis and back to G1)
Check telomere ttaggg longer than ccctaa for
the two G1 chromosomes (pter)
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k0=K.Chromosome(w, c)
k0.set_G1_to_G2()
k1 = k0.segregate()
#---------------------
ccctaa0 = k0.get_pter_G1(motif='ccctaa')
ttaggg0 = k0.get_pter_G1(motif='ttaggg')
#print ccctaa0, ttaggg0
ccctaa1 = k1.get_pter_G1(motif='ccctaa')
ttaggg1 = k1.get_pter_G1(motif='ttaggg')
#Check at pter
self.assertLess(ccctaa0,ttaggg0)
self.assertLess(ccctaa1,ttaggg1)
def test_4_telomere_G1_qter_GStrand(self):
""" Take a chromosome :
*G1 to G2
*segregate (mitosis)
Check telomere ttaggg longer than ccctaa(qter)
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k0=K.Chromosome(w, c)
k0.set_G1_to_G2()
k1 = k0.segregate()
#---------------------------------------
ccctaa0 = k0.get_qter_G1(motif='ccctaa')
ttaggg0 = k0.get_qter_G1(motif='ttaggg')
#print ccctaa0, ttaggg0
ccctaa1 = k1.get_qter_G1(motif='ccctaa')
ttaggg1 = k1.get_qter_G1(motif='ttaggg')
self.assertLess(ccctaa0,ttaggg0)
self.assertLess(ccctaa1,ttaggg1)
def test_5_telomere_G2_(self):
""" Take a chromosome :
*G1 to G2
* compute telomre lenghth on G2 chrom
Check that telomere length on both chromatids are similar
"""
w=K.Watson(ttaggg=10000, ccctaa = 5000)
c=K.Crick(ttaggg=5000, ccctaa = 10000)
k0=K.Chromosome(w, c)
k0.set_G1_to_G2()
motif='ccctaa'
pter = k0.get_pter_G2(motif)
p1 = pter[0]
p2 = pter[1]
self.assertAlmostEqual(p1,p2,delta=400)
qter = k0.get_qter_G2(motif)
q1 = qter[0]
q2 = qter[1]
self.assertAlmostEqual(q1,q2,delta=400)
motif='ttaggg'
pter = k0.get_pter_G2(motif)
p1 = pter[0]
p2 = pter[1]
self.assertAlmostEqual(p1,p2,delta=400)
qter = k0.get_qter_G2(motif)
q1 = qter[0]
q2 = qter[1]
self.assertAlmostEqual(q1,q2,delta=400)
def test_6_telomere_length_method(self):
'''Check that telomere_length() method gives the same results than
the other methods get_pter|qter_G1|G2(motif='ccctaa'|'ttaggg'),
with an easier to use interface.
'''
w = K.Watson(ttaggg=10000, ccctaa = 5000)
c = K.Crick(ttaggg=5000, ccctaa = 10000)
k0 = K.Chromosome(w, c)
L = k0.telomere_length(end='pter', motif='ccctaa')
self.assertEqual(L,10000)
L = k0.telomere_length(end='qter', motif='ccctaa')
self.assertEqual(L,5000)
L = k0.telomere_length(end='pter', motif='ttaggg')
self.assertEqual(L,10000)
L = k0.telomere_length(end='qter', motif='ttaggg')
self.assertEqual(L,5000)
#Switch to G2:Two chromatids Chromosome
k0.set_G1_to_G2()
self.assertIsInstance(k0.telomere_length(end='pter', motif='ccctaa'),tuple)
self.assertIsInstance(k0.telomere_length(end='qter', motif='ccctaa'),tuple)
p1, p2 = k0.telomere_length(end='pter', motif='ccctaa')
q1, q2 = k0.telomere_length(end='qter', motif='ccctaa')
self.assertAlmostEqual(p1,p2,delta=300)
self.assertAlmostEqual(q1,q2,delta=300)
#Switch from G2 to M:One chromatid Chromosome
k0.segregate()
self.assertIsInstance(k0.telomere_length(end='pter', motif='ccctaa'),int)
if __name__ == '__main__':
ut.main()
Raw
unittest_Genome.py
# -*- coding: utf-8 -*-
"""
Created on Thu Oct 11 14:20:03 2012
@author: Jean-Pat
"""
import unittest as ut
import ClassChromosome as K
import Genome as G
class Test_Genome(ut.TestCase):
""" class for unit tests for Genome class """
def test_genome01(self):
'''Check creation of a genome instance
'''
watson = K.Watson(ttaggg=1200, ccctaa=5500)
crick = K.Crick(ttaggg=5500, ccctaa=1200)
c1 = K.Chromosome(watson, crick, name='1', homolog='a')
c2 = K.Chromosome(homolog='b')
x = K.Chromosome(K.Watson(ttaggg=12000, ccctaa=4000),
K.Crick(ttaggg=4000, ccctaa=12000),
name='2', homolog='X')
y = K.Chromosome(K.Watson(ttaggg=20000, ccctaa=7000),
K.Crick(ttaggg=7000, ccctaa=20000),
name='2', homolog='Y')
g = G.Genome(c1, c2, x, y)
self.assertIsInstance(g, G.Genome)
def test_genome02(self):
''' Check can get the chromosomes name & homologous chrom name
'''
watson = K.Watson(ttaggg=1200, ccctaa=5500)
crick = K.Crick(ttaggg=5500, ccctaa=1200)
c1 = K.Chromosome(watson, crick, name='1', homolog='a')
c2 = K.Chromosome(homolog='b')
x = K.Chromosome(K.Watson(ttaggg=12000, ccctaa=4000),
K.Crick(ttaggg=4000, ccctaa=12000),
name='2', homolog='X')
y = K.Chromosome(K.Watson(ttaggg=20000, ccctaa=7000),
K.Crick(ttaggg=7000, ccctaa=20000),
name='2', homolog='Y')
g = G.Genome(c1, c2, x, y)
ploidy = g.get_chromosomes_list()
#print(ploidy)
self.assertEqual(list(ploidy), ['1','2'])
homologous = g.get_homologous_tags(pair='1')
self.assertEqual(homologous, ['a','b'])
homologous = g.get_homologous_tags(pair='2')
self.assertEqual(homologous, ['X','Y'])
def test_genome03(self):
'''Check can get chromosome phase
'''
watson=K.Watson(ttaggg = 1200, ccctaa = 5500)
crick=K.Crick(ttaggg = 5500, ccctaa = 1200)
c1 = K.Chromosome(watson, crick, name='1',homolog='a')
c2 = K.Chromosome(homolog='b')
g = G.Genome(c1, c2)
f1=g.get_chrom_phase(pair='1',hog='a')
self.assertEqual(f1,'G1')
f2=g.get_chrom_phase(pair='1',hog='b')
self.assertEqual(f1,'G1')
g.set_G1_to_G2()
f1=g.get_chrom_phase(pair='1',hog='a')
self.assertEqual(f1,'G2')
f2=g.get_chrom_phase(pair='1',hog='b')
self.assertEqual(f2,'G2')
def test_genome31(self):
'''check can get a pair of chrom with
'''
#TODO:Check that can get a pair of chromosomes with G.get...
pass
def test_genome04(self):
'''Check can get telomere length with self.get_telo_length method
'''
watson=K.Watson(ttaggg = 1200, ccctaa = 5500)
crick=K.Crick(ttaggg = 5500, ccctaa = 1200)
c1 = K.Chromosome(watson, crick, name='1',homolog='a')
c2 = K.Chromosome(homolog='b')
g = G.Genome(c1, c2)
p1a = g.get_telo_length()
p1a = g.get_telo_length(pair='1',hog='a', tel='pter',cycle='G1',motif='ccctaa')
self.assertEqual(p1a,1200)
def test_genome05(self):
'''Check can get telomere length with self.get_telomere_length method
'''
watson = K.Watson(ttaggg=1200, ccctaa=5500)
crick = K.Crick(ttaggg=5500, ccctaa=1200)
c1 = K.Chromosome(watson, crick, name='1', homolog='a')
c2 = K.Chromosome(homolog='b')
g = G.Genome(c1, c2)
p1a = g.get_telomere_length()
self.assertEqual(p1a, 1200)
q1a = g.get_telomere_length(tel='qter')
self.assertEqual(q1a, 5500)
def test_genome06(self):
'''Check telomere length difference between homologs
'''
watson = K.Watson(ttaggg=1200, ccctaa=5500)
crick = K.Crick(ttaggg=5500, ccctaa=1200)
wat1 = K.Watson(ttaggg=7500, ccctaa=9708)
crk1 = K.Crick(ttaggg=9708, ccctaa=7500)
c1 = K.Chromosome(watson, crick, name='1', homolog='a')
c2 = K.Chromosome(wat1, crk1, name='1', homolog='b')
g = G.Genome(c1, c2)
self.assertEqual(g.delta_telo(), 1200 - 7500)
self.assertEqual(g.delta_telo(tel='pter'), 1200 - 7500)
self.assertEqual(g.delta_telo(tel='qter'), 5500 - 9708)
def test_genome_state(self):
'''Check genome state:
'''
if __name__ == '__main__':
ut.main()
Raw
WatsonCrick.py
# -*- coding: utf-8 -*-
"""
Created on Wed Sep 5 11:05:37 2012
@author: Jean-Patrick Pommier
"""
import numpy as np
#The 5' ccctaa exonuclease deletion is simulated
#by a Binomial law(n=deletion, p=p_deletion)
#RNA = 20
#deletion = 400
#p_deletion = 0.40
#==============================================================================
# Watson & Crick seems to be the same class.but they are symetrical. A Watson
# object produces a Crick object and reciprocally as DNA replication does.
#==============================================================================
class Watson(object):
#Better to put these cte as class variables, same value for all instances
RNA = 20
deletion = 400
p_deletion = 0.40
def __init__(self, ttaggg = 5000, ccctaa = 10000):
self.__ttaggg = ttaggg
self.__ccctaa = ccctaa
self.in_use = False
def delete5prime(self):
#simulation action of 5' exonuclease (Makarov 1997)
self.__ccctaa = self.__ccctaa - (Watson.RNA+np.random.binomial(Watson.deletion,Watson.p_deletion))
def replicate(self):
self.__neottaggg = self.__ccctaa
self.__neoccctaa = self.__ttaggg
#make a complementary single strand
self.compstrand = Crick(self.__neottaggg,self.__neoccctaa)
#now delete Watson strand and Crick strand on 5' side
self.delete5prime()
self.compstrand.delete5prime()
return self.compstrand
def lengthen_ttaggg(self, l):
#use for telomerase
self.__ttaggg = self.__ttaggg+l
def lengthen_ccctaa(self, l):
#use for telomerase+DNA pol
self.__ccctaa = self.__ccctaa+l
def getTTAGGG(self):
return self.__ttaggg
def getCCCTAA(self):
return self.__ccctaa
ttaggg = property(fget=getTTAGGG, fset=None)
ccctaa = property(fget=getTTAGGG, fset=None)
class Crick(object):
RNA = 20
deletion = 400
p_deletion = 0.40
def __init__(self,ttaggg = 10000, ccctaa = 5000):
self.__ttaggg = ttaggg
self.__ccctaa = ccctaa
self.in_use = False
def delete5prime(self):
self.__ccctaa = self.__ccctaa - (Crick.RNA+np.random.binomial(Crick.deletion,Crick.p_deletion))
def replicate(self):
self.__neottaggg = self.__ccctaa
self.__neoccctaa = self.__ttaggg
#make a complementary single strand
self.compstrand = Watson(self.__neottaggg,self.__neoccctaa)
#now delete Watson strand and Crick strand on 5' side
self.delete5prime()
self.compstrand.delete5prime()
return self.compstrand
def lengthen_ttaggg(self, l):
#use for telomerase
self.__ttaggg = self.__ttaggg+l
def lengthen_ccctaa(self, l):
#use for telomerase+DNA pol
self.__ccctaa = self.__ccctaa+l
def getTTAGGG(self):
return self.__ttaggg
def getCCCTAA(self):
return self.__ccctaa
ttaggg = property(fget = getTTAGGG, fset = None)
ccctaa = property(fget = getTTAGGG, fset = None)
#==============================================================================
#
#==============================================================================
def test_watson_replication():
w = Watson()
print( 'w TTAGGG:',w.getTTAGGG())
print( 'w CCCTAA:',w.getCCCTAA())
neoc = w.replicate()
print( 'w post replication CCCTAA:',w.getCCCTAA())
print( 'Crick strand post rep, CCCTAA:',neoc.getCCCTAA())
w.__ccctaa = 200
print( w.getCCCTAA())
if __name__ == "__main__":
test_watson_replication()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment