Modeling the Telomere length distribution with Chromosomal Segregation along Synchronous Cellular Divisions

Telomeres: structure and functions

If it is necessary to explain what the telomeres are, let's watch some lectures:

First from Elisabeth Blackburn:

 

Or from Titia de Lange:

 

Modeling Telomere Dynamic:

Levy et al. in Telomere End-replication Problem and Cell Aging, modeled the distribution of telomeric erosion with cells divisions at a single telomere. The model accounts for the decrease in somatic cell division capabilities with in-vitro passages of fibroblasts (Hayflick's limit).

Here, the telomere length of several chromosomes is modeled numerically with segregation of chromosomes in daughter cells. The model takes the 5' degradation of the CCCTAA strand into account. The aim of the simulation is to compare the distributions with QFISH data. The model shows that the telomeres length at pter and qter are correlated for a given homologous chromosome. It can be used to derive statistical tests to detect telomere length difference between homologous chromosomes.

Code Installation

Download several python modules and the jupyter notebook in the same directory:

• WatsonCrick.py
• ClassChromosome.py
• Genome.py
• Cellular.py
• telomere-length-distribution-in-synchronous-dividing-cells-ipynb

Some unit tests can be run, a jupyter notebook is available (see the link to gist), they yield:

Model structure:

 The  model consists in a class Cell which have a Genome, which have chromosomes, which have two complementary single strands DNA (Watson/Crick). Each strand has two telomeres.

A single strand object can be instantiated as follow:

A chromosome object is instantiated with two complementary Watson/Crick strands:

For example, two chromosomes, rosalind and franklin can be instantiated as follow:

A single chromatid chromosome, in the G1 state of the cell cycle, can be triggered to the G2 state:

After a round of DNA replication, 5' CCCTAA motifs are incompletely replicated and TTAGGG 3' can be randomly degraded, leading to shorter telomeres on metaphasic chromosomes on both chromatids at pter and qter. Iin the G2 state a chromosome has two chromatids, the length of a telomere on a given arm (pter ou qter) is given by two values:

The mitosis is triggered by chromosome segregation:

The ros chromosome object keeps its signature before and after the mitosis:

Signature of the ros chromosome in the G2 state before mitosis

Definition of a Genome object:

It is convenient to define a genome object, it is instantiated for example with two chromosomes (2N=2):

The genome state (G1/G2), the chromosomes of the genome can be accessed to read the telomere length (CCCTAA or TTAGGG motifs), the length difference between homologous telomeres:

Cell object has a genome:

Unit test results for Cell object

Once a single genome object is instantiated, a single cell object can be instantiated. Let's build a cell with 2N=10:

The length of the telomeres (bp) in that single cell can be read from a pandas dataframe:

Let's allow 10 synchronous cellular divisions from 1 cell to 1024 cells. Then as for fibroblasts cultures, let's make passages, that is one allows only half randomly chosen cells to divide once:

Results

The envelope of the telomere length distribution at a single telomere after 14 divisions (10 divisions + 4 passages) seems to be gaussian:

Correlation between the length at pter and qter on the same homolog

The length of the telomeres (pter, qter) of the same homologous chromosome (the maternal one for example) are correlated:

It is possible to plot the telomeres length of two homologous chromosomes:

The length of the telomeres belonging to different homologous are not correlated.

Let's plot the length of the telomere from the paternal chromosome 1 at qter as a function of the length of the telomere at 1qter on the maternal homolog:

Mean telomere length in synchronous dividing cells population:

In telomerase negative cells as modeled here, the telomere length decreases with cells divisions.

The decrease of telomere length depends on the 5' exonuclease activity. The amount of degraded DNA is modeled here by random variable, X, following a binomial law (N=400, p=0.4) :

Distribution of X: P(130<X<200)=0.999

to fit the data of Makarov et al. (1997) as follow:

 The expected value of X, the deletion length, is :

taking the RNA primer deletion (20 bp), the mean length of the 3' s G strand is:

160 +20 bp

As published in 1992, the standard deviation of the telomere length increase at each cell division.

So the mean telomere length, $mean+/-std$, was plotted as a function of the cell divisions after an expansion of ten cellular divisions and four passages.

The mean telomere length decreases by 90 bp/div

and the heterogeneity increases:

At each passage the length of each telomere (column) of each cell (row) of the simulation is copied in a pandas data-frame, for example at passage P4:

 Then the mean telomere length in the whole cells population can be calculated:

Using seaborn, the telomere length of some homologs can be compared:

Population of Mixed Cells

Cells from different passages can be mixed, for example in equal proportions. The correlation between pter and qter of a given homolog should start to vanish:

Senescence in-silico:

Short telomeres trigger an irreversible transition to $G_0$ state of the cell cycle (senescence). Let's take a cells population, make successive passages. Initial cells population is expended from one cell (telomerase off) after eight  divisions. The genome of the initial cell is instantiated with 2N=4 and with the length of the shortest telomere set to 2000 bp.

 

At each passage, it's possible to count the cells in the $G_0$ state or to calculate the confluence of the cells population. The threshold telomere length triggering the G0 state is set arbitrarily to 200 bp.  The initial population is passed 20 times:

The confluence is calculated as the number of cells at a given passage by the initial number of cells. With the shortest initial telomere length set to 2000 bp, the confluence reach 50% after 18 doublings:

Senescence in silico: confluence  (100 N cell/ N initial cells) decreases with cells divisions

Python modules and jupyter notebook

	# -*- 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

view raw Cellular.py hosted with ❤ by GitHub

	# -*- 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()

view raw ClassChromosome.py hosted with ❤ by GitHub

	# -*- 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)

view raw Genome.py hosted with ❤ by GitHub

Loading

Sorry, something went wrong. Reload?

Sorry, we cannot display this file.

Sorry, this file is invalid so it cannot be displayed.

view raw Telomere Length Distribution in Synchronous Dividing Cells.ipynb hosted with ❤ by GitHub

	# -*- 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()

view raw unittest-watsoncrick.py hosted with ❤ by GitHub

	# -*- 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()

view raw unittest_Cell.py hosted with ❤ by GitHub

	# -*- 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()

view raw unittest_Chromosome.py hosted with ❤ by GitHub

	# -*- 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()

view raw unittest_Genome.py hosted with ❤ by GitHub

	# -*- 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()

view raw WatsonCrick.py hosted with ❤ by GitHub