| |
CyberFly: Genetics Simulatior
by Lex Berezhny and David Massant
Introduction
The Drosophila melanogaster fruit fly is a little insect about 3mm long, that typically accumulates around spoiled fruit. It is also one of the most valuable organisms in biological research, particularly in genetics and developmental biology. Drosophila has been used as a model organism for research for almost a century, and today, several thousand scientists are working on many different aspects of the fruit fly.
Part of the reason people work on it is historical - so much is already known about it that it is easy to handle and well-understood - and part of it is practical its a small animal, with a short life cycle of just two weeks, and is cheap and easy to keep large numbers. Mutant flies, with defects in any of several thousand genes, are available, and the entire genome has recently been sequenced. For this reason, we were introduced to the Drosophila in high school labs where we were taught genetics by repeating a famous experiment done in the early twentieth century by Thomas Hunt Morgan. This experiment also introduced us to the negative aspects involved in working with the flies. The anesthetic used on the flies began to have a similar affect on the class (augmented by the tedium of counting the flies one by one under a microscope), which is dangerous given that inattention could result in a pair of flies escaping and causing a huge infestation in the labs.
Morgan's experiment, which so many of our classmates enjoyed, was the first to link a gene to a particular chromosome. This was the discovery of the mechanism for sex linkage. In his experiment, Morgan mated a group of wild type (those which occur in the wild; no mutations; red eyed) females with a group of white eyed males. In fruit flies, white eyes is a mutation (all fruit fly mutations are recessive) that we now know only occurs on the X chromosome. This means that the two X chromosomes in the females had the wild type gene (w+) and the one X chromosome in the males had the mutation (w). In the first generation he found all the flies to have red eyes (this meant that all the females were heterozygous for the mutation and all the males had an X chromosome without the mutation). After mating this next generation, Morgan found no female flies with the white eye mutation and half the population of male flies with the white eye mutation. Thus half the males had gotten an X chromosome with the mutations and the other half had gotten an X chromosome without the mutation. This experiment conclusively showed that the white eye gene occurs only on the X chromosome.
A great deal of work in testing for autosomal linkage in the mating of fruit flies has been done as well. Autosomal linkage is used to describe genes that are located on the same chromosome and thus travel together. An example in humans is the gene that codes for freckles and the gene that codes for red hair. While not all people with red hair have freckles (or vice versa), the majority of people with red hair also have freckles. Fruit flies have four chromosomes (compared to a humans twenty-six) each of which has a host of non-lethal mutations.
There is another biological process that occurs during gamete formation that results in autosomal linkages being broken. This process, called crossing over, is the swapping of genetic material between the like chromosomes of an individual. The chances of a cross over can be calculated using the actual physical distances between the autosomally linked genes.
These genetic processes have the overall effect of complicating any effort to calculate the expected progeny of any given fruit flies. The objective of our project is to avoid the hassle of dealing with actual flies and bypass complex math by creating a computer program to make the flies for us. This program would be able to emulate the progeny of fruit flies given the characteristics and genotype of the parents accounting for autosomal linkage, sex linkage, and crossing over.
The software to address this problem is called CyberFly and was written in Python. Python is a high level, dynamically typed, interpreted, programming language. Similar to its counterparts like the industry standard languages Java and C++, it has full object-oriented programming support providing a facility for creating efficient high-level data structures, encapsulated behavior and black-box implementations.
Dynamic typing and its interpreted nature set Python apart from many other languages. Neither Java nor C++ offer such features. Dynamic typing allows the programmers to stay concentrated on the problem being addressed instead of meticulously worrying about the details of how their data will be represented inside the computer. Dynamic typing will do this automatically.
Unlike C++ and Java, which are compiled languages, Python is interpreted. This removes the extra step of having to recompile code upon every change. Further, the interactive interpreter, available due to Pythons interpreted nature, allows for quick experimentation of specific algorithms without the overhead of weaving it into the software. Many algorithms can be tested quickly before the best one is chosen for inclusion into the actual software. The benefits of reduced development time and increased programmer productivity almost always outweigh the runtime speed constraint of interpreted languages. CyberFly is an example where speed and efficiency are second to ease of development and development time constraints.
Python was conceived by Guido van Rossum in 1990 at CWI in Amsterdam initially for educational and scientific purposes. The language has continued to grow and foster in those communities, gaining much popularity at universities and research centers. The National Air and Space Administration is just one organization that uses Python extensively for work flow and personnel management software, testing shuttle equipment, integrating heterogenous systems and developing mission critical software. Recently, Pythons popularity has spilled over into to the commercial arena, extending the already large Python community and adding to Pythons rich set of tools and software. Industrial Light and Magic, the special effects company responsible for many movies including Star Wars, uses Python to prototype sophisticated modeling software and to integrate legacy applications with newer or in-house developed systems.
Procedures
CyberFly was born on the back of a coffee stained napkin. A symbol that never left the project as the engineering process unfolded consuming most of the alloted project time and exceeding double the programming time. Much of the initial logic and algorithms had been designed on paper prior to writing a single line of code.
Sketching a high level object model was the initial stage in setting the complexity and problem domain of the genetics and simulation engine in CyberFly. Once the birds eye view was taken, certain aspect where magnified for clarity, and as more implementation details where discussed and tested via Pythons interactive interpreter the object model changed shape to suit the new requirements and observations.
Once enough notes, diagrams and sample code where produced the focus moved from the design stage to synthesizing code. Initially the coding process was trivial as most of the logic and algorithms had been planned out. Connecting the different pieces, algorithms and functions was less trivial but using Pythons interactive interpreter it was possible to mix and match parts and find elegant solutions to defining the flow of execution.
Conceptually and practically the project was divided into two portions. The genelib library contains all required genetics and simulation engine logic. On top of genelib sits the GUI or graphical user interface layer which was written using the wxWindows graphics toolkit.
The genelib library consists of four top level objects, namely Gene, Chromosome, GeneMap and Fly. The atomic Gene object provides a contained 'blackbox' interface to the details of a gene, including it's name, chromsome and chromosomal position. Gene objects are commonly contained within Chromosome objects. Chromosome objects are more advanced in that they can be asked to cross-over with another chromosome. Finally a Fly object is responsible for maintaining a collection of chromosomes and recursively crossing them when asked to mate.
The genelib was designed with the intent of having an easy to use API. There is a trivial number of steps to run the simulation. An initial setup process includes loading the drosophila gene map file into the GeneMap object, and then instantiating a male and female Fly objects. The pair of Fly objects are then assigned the requested genotypes and become ready for mating.
Every instance of the Fly object has an intuitive mate method which is executed by passing two arguments, respectively: the partner fly object to mate with and the number of baby flies to generate. The mate method then loops performing parental crossing-over to make each baby, accounting for autosomal and sex linkage.
Examining this in detail, top down, the mate method actually acts as a trigger for other auxiliary methods and is only responsible for calling those methods and keeping track of the generated babies. The cross_n_choose method, spawned by mate, works by looping over pairs of paternal and maternal chromosomes and calling the cross method of one of the chromosomes while passing the other as an argument. The result of this call returns two crossed chromosomes, one of which is randomly chosen for the baby.
The complexity of the cross method is daunting but serves the essence of this entire project. During the initially design process several goals, even limitations, where set to guarantee a timely completed product. The biggest issue was deciding how general the crossing-over algorithms should be, complexity went up as large portions of the algorithms were generalized and eventually several limitations where decided upon. The main limitation is to only allow single and double crossing, another is an inability to deal with more than four traits.
The types of crossing-over calculations performed are based on the number of genes on the chromosome. For two and three genes on the chromosome only the single cross is calculated and if there are four genes both the single and double cross calculations are performed. Two private functions of the Chromosome class _single_cross and _double_cross are used for the calculations (see Appendix A). Based on the calculations the cross method moves the genes around emulating real life crossing-over and finally returns two inter-crossed chromosomes.
Once the list of Fly objects is returned from a mate method call it is up to the GUI to render the results and perform statistical analysis on the raw data.
Results
CyberFly uses Pythons built-in random module for adding variability to the crossing-over process to further enhance the simulation. Python uses the standard Wichmann-Hill generator, combining three pure multiplicative congruential generators of modulus 30269, 30307 and 30323 (Wichmann). Its period of 6,953,607,871,644 is more than ideal for our purposes. On average there are 24 calls into the random module per baby, thus generating ten thousand babies during a simulation would mean ~240,000 calls. The simulation could be run 28,973,366 with ten thousand babies each time before one might see duplicate results.
Conclusion
There are many extensions to this project that would make it better suited for more complicated simulations. One option that we include is the changing of the genemap. In our program we use the genemap of the fruit fly (see Appendix B), however it is possible to load and soon to create new genemaps. This would allow our program to run using the genetic make up of any animal.
Another possibility for our program would be to allow it to accept data and build a genemap from that data using similar statistics. This option would probably have the most potential especially for researchers, however students would also learn a great deal from creating their own genemaps. If the exact sequence data were entered in place of the genemap, the program could actually calculate the autosomal linkage across base pairs (much shorter than the distances between genes).
Runtime speed has not been a critical problem but there is potential for serious speed degradation if more features and processor intensive calculations are added. Fortunately Python is easily extensible from the C programming language and there is room for reimplementing many of the algorithms and processor intensive tasks using C. The speed increase could be very substantial. C is a very low level language that is one step beyond assembly programming, in fact Python was written in C. Thus it is very likely that some portions of CyberFly will be rewritten in C for speed increase.
Bibliography
Wichmann, B. A. & Hill, I. D., ``Algorithm AS 183: An efficient and portable pseudo-random number generator'', Applied Statistics 31 (1982) 188-190.
Appendix A
import string
from random import random, randrange
DEBUG = 0
class Gene:
"""
Represents a single gene. Gene objects
are usually children of the Chromosome
object.
Special properties of a gene object:
subtracting:
you can subtract two genes to get
the distance between them.
sorting:
you can sort a list of genes which
will put them in positional order
"""
def __init__(self, chrom, gene, pos):
""" Represents a single gene. """
self.chrom = chrom # name of chromosome (string)
self.gene = gene # name of gene (string)
self.pos = pos # position (float)
### Get Methods ##############################
def getChrom(self):
return self.chrom
def getName(self):
return self.gene
def getPosition(self):
return self.pos
### Overloaded Basic Opeartors ###############
def __str__(self):
return self.getName()
__repr__ = __str__
def __sub__(self, other):
diff = abs(float(self.pos) - float(other.pos))
if diff > 50:
diff = 50
return diff
def __int__(self):
return int(self.pos)
def __float__(self):
return float(self.pos)
### Overload Compare/Sorting Operator ########
def __cmp__(self, other):
""" Used to sort a list of genes in increasing
positional order.
"""
if self.pos == other.pos:
return 0
if self.pos < other.pos:
return -1
if self.pos > other.pos:
return 1
class Chromosome:
"""
A chromosome object contains a list of genes
that comprise it and is able to perform
crossing-over with another chromosomes genes.
"""
def __init__(self, name, genes = []):
self.name = name
self.genes = genes[:]
self.probs = [] #probabilities
### Set Methods ##############################
def append(self, gene):
self.genes.append(gene)
### Get Methods ##############################
def phenotype(self, other):
""" Compares two chromosomes looking
for identical genes.
"""
homo = []
for gene in self.genes:
if gene in other.genes:
homo.append(gene)
return homo
def __repr__(self):
self.genes.sort()
return "\n"+self.name+"\n"+string.join(map(str, self.genes), "\n")
### Crossing-over Methods ####################
def _chance(self, base, other):
""" Chance of an event occuring. """
prob = 1
debug_eqs = []
for gene in other:
prob *= (gene-base)/100 #devides the difference between two genes by 100
debug_eqs.append("((%.1f-%.1f)/100)" % (gene, base))
if DEBUG > 4: print string.join(debug_eqs, "*"),"=",prob
return prob*100
def _single_cross(self):
""" Compiles the chances of single crosses. """
sum = 0
probs = []
for i in range(0, len(self.genes)):
sum += self._chance(self.genes[i], self.genes[:i]+self.genes[i+1:])
probs.append(sum)
return probs
def _double_cross(self, offset):
""" Chances of doing a double cross.
"""
sums = []
res = ((self.genes[0]-self.genes[2])*\
(self.genes[1]-self.genes[2])*\
(self.genes[0]-self.genes[3])*\
(self.genes[1]-self.genes[3]))/(100**3)
sum = res+offset
sums.append(sum)
res = ((self.genes[0]-self.genes[1])*\
(self.genes[1]-self.genes[3])*\
(self.genes[0]-self.genes[2])*\
(self.genes[2]-self.genes[3]))/(100**3)
sum += res
sums.append(sum)
res = ((self.genes[0]-self.genes[1])*\
(self.genes[0]-self.genes[3])*\
(self.genes[1]-self.genes[2])*\
(self.genes[2]-self.genes[3]))/(100**3)
sum += res
sums.append(sum)
return sums
def _calc_cross(self):
""" Calculate the crossing-over probabilities.
"""
gene_count = len(self.genes)
if gene_count == 2 or gene_count == 3:
return self._single_cross()
elif gene_count == 4:
single = self._single_cross()
double = self._double_cross(single[-1])
single.extend(double)
return single
return [] #no crossing over occurs
def cross(self, pair):
""" Crosses this chromosome with
it's pair.
"""
# make sure the first one is ordered
# by increasing position
# then make a deep copy of both
self.genes.sort()
set1 = self.genes[:]
set2 = pair.genes[:]
probs = self._calc_cross()
rand = random()*100
#switch genes
def switch(a,b,gene):
""" Cases:
If in both, do nothing.
If in a, remove from a and append to b.
"""
if not gene in b:
a.remove(gene)
b.append(gene)
# decide case
case = -1
for i in range(len(probs)):
if probs[i] >= rand:
case = i
break
if case == 0:
if DEBUG > 2: print "Case I"
switch(set1, set2, set1[0])
elif case == 1:
if DEBUG > 2: print "Case II"
switch(set1, set2, set1[1])
elif case == 2:
if DEBUG > 2: print "Case III"
switch(set1, set2, set1[2])
elif case == 3:
if DEBUG > 2: print "Case IV"
switch(set1, set2, set1[3])
elif case == 4:
if DEBUG > 2: print "Case V"
one, two = set1[2], set1[3]
switch(set1, set2, one)
switch(set1, set2, two)
elif case == 5:
if DEBUG > 2: print "Case VI"
one, two = set1[0], set1[3]
switch(set1, set2, one)
switch(set1, set2, two)
elif case == 6:
if DEBUG > 2: print "Case VII"
one, two = set1[1], set1[3]
switch(set1, set2, one)
switch(set1, set2, two)
else:
pass
return Chromosome(self.name,set1), Chromosome(self.name,set2)
class GeneMap:
""" A pool of genes. """
def __init__(self, map_file_name):
""" Parses chromosome definition file into memory. """
self.genes = []
self.chromes = []
chromosome = ""
for line in open(map_file_name).readlines():
line = line.strip()
if not line:
chromosome = ""
else:
if not chromosome:
chromosome = line
self.chromes.append(chromosome)
else:
# TODO: handle exceptions for malformed data (ie, distance must be float)
distance, gene = line.split('\t')[0], string.join(line.split('\t')[1:])
self.genes.append(Gene(chromosome, gene, float(distance)))
### Get Methods ##############################
def getGene(self, name):
""" Tries to find a gene object with
corresponding name.
Returns None otherwise.
"""
for gene in self.genes:
if name == gene.getName():
return gene
return None
def getChroms(self):
return self.chromes
def getGenes(self, type):
lst = []
for g in self.genes:
if g.getName().find(type) > -1:
lst.append(g)
return lst
def __repr__(self):
for gene in self.genes:
print gene
class Fly:
def __init__(self, map):
""" Instantiates fly object. """
self.map = map
# chromosome sets
self.set1 = {}
self.set2 = {}
for chrom in self.map.getChroms():
self.set1[chrom] = Chromosome(chrom)
self.set2[chrom] = Chromosome(chrom)
def addGene(self, name, genotype):
""" Adds gene to fly.
name - unique name of gene
genotype -
hw (homozygous wild) ignored
ht (heterozygous) set1
hm (homozygous normal) set1,set2
"""
gene = self.map.getGene(name)
if genotype == "hw":
pass
elif genotype == "ht":
self.set1[gene.getChrom()].append(gene)
elif genotype == "hm":
self.set1[gene.getChrom()].append(gene)
if self.sex != 'male' and gene.getChrom()[0] != 'X':
self.set2[gene.getChrom()].append(gene)
def cross_n_choose(self):
""" Crosses it's own chromosomes and returns
one of the two randomly selected
chromosomes.
"""
chroms = self.set1.keys()
chroms.sort()
set = {}
for chrom in chroms:
c = self.set1[chrom].cross(self.set2[chrom])[randrange(0,2)]
set[chrom] = c
return set
def mate(self, other, n_babies):
""" Mate this fruitfly with the 'other' one. """
babies = []
while n_babies:
baby = Fly(self.map)
baby.set1 = self.cross_n_choose()
baby.set2 = other.cross_n_choose()
baby.sex = ['female', 'male'][randrange(0,2)]
babies.append(baby)
n_babies -= 1
return babies
def phenotype(self, use_cached = 1):
""" Returns the phenotype. """
if use_cached and hasattr(self, "pheno"):
return self.pheno
chroms = self.set1.keys()
chroms.sort()
self.pheno = []
for chrom in chroms:
self.pheno.extend(self.set1[chrom].phenotype(self.set2[chrom]))
self.pheno.sort()
return self.pheno
def pprint(self):
""" Generate phenotype for fly. """
chroms = self.set1.keys()
for x in chroms:
print self.set1[x]
Appendix B
X chromosome 1
0.0 yellow body
0.0 scute bristles
1.5 white eyes
3.0 facet eyes
5.5 echinus eyes
7.5 ruby eyes
13.7 crossveinless wings
20.0 cut wings
21.0 singed bristles
27.7 lozenge eyes
33.0 vermillion eyes
36.1 miniature wings
43.0 sable body
44.0 garnet eyes
56.7 forked bristles
57.0 bar eyes
62.5 carnation eyes
body chromosome 2
1.3 star eyes
4.0 held-out wings
13.0 dumpy wings
16.5 clot eyes
48.5 black body
51.0 reduced bristles
54.5 purple eyes
54.8 short bristles
55.0 light eyes
57.5 cinnabar eyes
66.7 cabrous eyes
67.0 vestigial wings
72.0 lobe eyes
75.5 curved wings
100.5 plexus wings
104.5 brown eyes
107.0 speck body
body chromosome 3
0.0 roughoid eyes
19.2 javelin bristles
26.0 sepia eyes
26.5 hairy body
41.0 dichaete bristles
44.0 scarlet eyes
48.0 pink eyes
50.0 curled wings
58.2 stubble bristles
58.5 spineless bristles
58.7 bithorax body
62.0 stripe body
63.0 glass eyes
69.5 hairless bristles
70.7 ebony body
74.7 cardinal eyes
91.1 rough eyes
100.7 claret eyes
106.2 minute bristles
body chromosome 4
0.0 bent wing
0.0 eyeless
|
|