Sunday, June 19, 2016

The Pursuit of Happiness or Happyness

Most people want to be happy. Except for the rare idiosyncratic individual the world really is a place filled with people who are in the pursuit of happiness. Happiness is a weird and sometimes illusive thing. The movie about The Pursuit of Happiness does a good job of capturing the illusive nature of Happyness.

Would eating the Gumbo make you happy? It has a history as much as Thomas Jefferson does, or may be even more.
How about modelling the evolution of directed graphs along a phylogeny? Would it make any difference what those directed graphs represent? They could be anything really, a map of the brain (connectome), parasite transmission networks, protein protein interaction networks or simply a simulation.  
In the end we are still only chasing after it.

Who has found it?

The unlikeliest of people have found it.  No it has nothing to do with Pittsburg skyline and everything to do with the mind.

Monday, May 30, 2016

Dolphins, crows and apes - as clever as it gets

Dolphins, just like crows & apes are very smart. In some way they represent the pinnacle of the independent evolution of intelligence in species that dwell in the oceans, fly in the sky and walk on land respectively. What can the brains of these distinct yet similar taxa tell us about intelligence? Will they be able to provide crucial insight needed to understand intelligence, thought and the brain? How can they guide artificial intelligence research?

All these may seem far fetched questions for another day. However, we are not too far away. A study published in the year 2013, "Large-scale network organization in the avian forebrain: a connectivity matrix and theoretical analysis" was able to generate a preliminary map of a avian forebrain. Work in this field is progressing at an incredible pace. It might be worth noting that one of the co-authors Murray Shanahan is actually a Professor of Cognitive Robotics and has written the book "The Technological Singularity". So the fields are not so far apart after all. 

Few months ago, while writing up my PhD thesis "Speciation genomics: A perspective from vertebrate systems" i began to realize how intricately linked the world is. One needs to understand the evolutionary genetics of phenotypic traits to be able understand speciation and adaptation. This understanding of genetics will play an important role someday in the future to look at traits like "intelligence". We may infact be able to unravel the great mysteries of the brain and its evolution. 

Understanding the brain, its evolution and genetics definitely have their own merits. The next step into the world of artificial intelligence and culture actually seems exciting at this point. What kind of morality would different machine cultures create? Would speciation "co-evolve" with culture in artificial systems as much as it seems to in the natural world?

Sunday, May 1, 2016

Ensembl Perl API to get all intron lengths in Human genom

The below script will get all stable id's from Ensembl and prints out the intron lengths for each transcript of every gene. Along with intron length, the flanking exon id's are also printed. One can get the upstream and downstream intron length for each exon using the output.

The output would look like this:
Gene Id                     Transcript Id            Previous Exon          Next Exon               Intron length
ENSG00000084674 ENST00000233242 ENSE00000932268 ENSE00000932269 717 ENSG00000084674 ENST00000233242 ENSE00000932269 ENSE00000932270 2338 ENSG00000084674 ENST00000233242 ENSE00000932270 ENSE00000932271 112 ENSG00000084674 ENST00000233242 ENSE00000932271 ENSE00000719046 1100 ENSG00000084674 ENST00000233242 ENSE00000719046 ENSE00000718984 261 ENSG00000084674 ENST00000233242 ENSE00000718984 ENSE00000932272 863 ENSG00000084674 ENST00000233242 ENSE00000932272 ENSE00000932273 1663 ENSG00000084674 ENST00000233242 ENSE00000932273 ENSE00000718481 1240 ENSG00000084674 ENST00000233242 ENSE00000718481 ENSE00000542194 482

 #!/usr/bin/perl  
 use strict;  
 use warnings;  
 use Bio::EnsEMBL::Registry;  
 use Bio::SeqIO;  
 use Getopt::Long;  
 my $registry = 'Bio::EnsEMBL::Registry';  
 ## Load the databases into the registry  
 $registry->load_registry_from_db(  
  -host => 'ensembldb.ensembl.org',  
  -user => 'anonymous'  
 );  
 ## Get the gene adaptor for human  
     my $gene_adaptor = $registry->get_adaptor( 'Human', 'Core', 'Gene' );  
     # Fetch my gene of interest usning ensemble ID  
     my @gene_ids = @{$gene_adaptor->list_stable_ids()};  
 foreach my $geneid(@gene_ids){  
 #print "$geneid\n";  
 my $gene = $gene_adaptor->fetch_by_stable_id($geneid);  
  foreach my $transcript (@{ $gene->get_all_Transcripts }) {  
   foreach my $intron (@{ $transcript->get_all_Introns }) {  
   print $gene->stable_id,"\t",$transcript->stable_id,"\t",$intron->prev_Exon->stable_id,"\t",$intron->next_Exon->stable_id,"\t",$intron->length,"\n";  
   }  
  }  
 }  

Wednesday, April 27, 2016

Visualizing isoforms of a gene as a undirected graph using sna package

The sna (social network analysis) package in R provides an easy to use interface for handling network data structures. Apart from numerous statistics that can be calculated on the graph, it is possible to visualize the graph using the "gplot" command.

Here, we use a perl script to convert the information regarding exon positions along the gene and transcript structure into the nos format. We are then able to see that the gene glutamate-cysteine ligase, catalytic subunit (GCLC) is actually made up of 5 different components.

Perl script to write the gene in nos format: (GeneRanked_exons contains list of exons with their positional rank in the gene. test.exon.order contains the list of exons in each transcript ordered by the exon positional rank in the transcript)
 #!/usr/bin/perl  
 #  
 ##perl printNetwork.pl GeneRanked_exons.txt test.exon.order > ENSG00000001084_graph.txt  
 #  
 my %exons=();  
 open(FILE, $ARGV[0]);  
 while (my $row = <FILE>) {  
  chomp $row;  
  @values=split(' ',$row);  
  $exons{$values[0]}=$values[2];  
 }  
 close FILE;  
 my %matrix;  
 my $maxexon=0;  
 open(FILE2, $ARGV[1]);  
 my $row = <FILE2>;  
 chomp $row;@values=split('\t',$row);$exoncount1=$exons{$values[0]};$trans1=$values[2];$trancount1=$values[3];  
 while (my $row = <FILE2>) {  
  if($exoncount1>$maxexon){$maxexon=$exoncount1;}  
  chomp $row;@values=split('\t',$row);$exoncount2=$exons{$values[0]};$trans2=$values[2];$trancount2=$values[3];  
  #print "$exoncount1\t$trans1\t$exoncount2\t$trans2\n";  
  if($exoncount2>$maxexon){$maxexon=$exoncount2;}  
  $ftrancount2=$trancount2-1;  
  #print "$exoncount1\t$trans1\t$exoncount2\t$trans2\t$trancount1\t$trancount2\n";  
  if(($trans1 =~ m/^$trans2$/i)&&($trancount1==$ftrancount2)){ $matrix{$exoncount1}{$exoncount2}=1;}  
  $exoncount1=$exoncount2;$trans1=$trans2;$trancount1=$trancount2;  
 }  
 print "1\n";  
 print "$maxexon $maxexon\n";  
 for ($i=0;$i<$maxexon;$i++){  
 print "0 0 ";  
 }  
 print "\n";  
 for ($k=1;$k<=$maxexon;$k++){  
      for ($i=1;$i<=$maxexon;$i++){  
      if(exists $matrix{$k}{$i}){$j=$matrix{$k}{$i};}  
  else{$j=0;}  
      print "$j ";  
      }  
      print "\n";  
 }  
 #     #1  
 #     #4 4  
 #     #0 0 0 0 0 0 0 0  
 #     #0 1 0 0  
 #     #0 0 1 1  
 #     #0 1 0 0  
 #     #0 0 1 0  
 #  
 #  
The output from the printNetwork command can be read into R using the read.nos function of the sna package. The graph can then be visualized using below lines in R.

 file1<-"ENSG00000001084_graph.txt"  
 library(sna)  
 read.nos(file=file1)->g  
 jpeg(paste(file1,".jpeg",sep=""))  
 gplot(g,gmode="graph",displaylabels=TRUE,main=file1)  
 dev.off()  
This produces a graph that looks like this: (each red dot is an exon with the number beside it being its positional rank).
While the actual gene on ensemble looks like this:


Saturday, December 5, 2015

Hawaizaada - could the first human controlled powered airplane have taken flight on the beaches of Mumbai?

A recent movie Hawaizaada tries to think of an alternative history. One in which the very first human controlled powered airplane takes flight not on the beaches near Kitty Hawk but the one in Mumbai, India. Would it have been possible for such a character to even exist in pre-independence India? Most definitely, the example of CV Raman comes to mind. Was the environment conducive for innovation? Was it the startup capital :) The very presence of a few stars gives hope for the future.

The ridicule faced by the main character when he tries to convince the Maharaja of Baroda is probably not very different from what was experienced by the people who tried to claim that the book "Vaimānika Shāstra" was a scientific treatise. While it is now established by the leading scientists of today that the book is no more than a bunch of fantasies that have no grounding in reality, one has to at the least appreciate that somebody thought of concocting such brilliant lies. The movie Hawaizaada tries to give color to the fantasies found in the book and follows the life of Subbaraya Shastry who tried to make these models with help from Dr.Talpade (who has been replaced by a school dropout in the movie) of Bombay. In reality none of these models ever worked, obviously due to their mythical nature. 

While the movie is a commendable attempt at trying to create a "new" genre in Indian Cinema, it was a box office failure which received mixed criticism, most of it being negative. This seems more of a comment on the state & taste of the Indian entertainment industry than on the quality of the movie. One can only hope that with increasing literacy and as understanding of the Jaredness of the world increases, such dreamers will flourish. 

Monday, November 23, 2015

Have empires become larger over time?

The list of largest Empires at Wikipedia tries to provide an unbiased list of empires, their size and the approximate time during which these empires existed. Even though the very definition of the term "empire" is fraught with controversies, an operational definition that one might accept would involve large tracts of land, resources and people that have been brought under the control of some sort of central political entity. 

Obviously, over time the world has become more connected and has required larger regions to be under a single dominion to be considered an empire. In the distant past, the control of a collection of tribes might have been enough to be considered an empire. Overtime, this has grown to require the control of more than one continent. Thus, in a way the very way that empires have been defined, requires that they become larger over time. So, one would expect a positive correlation between the size of the empires and the time period in which they have existed. 

We use the list compiled by wikipedia from various scholarly sources to test this hypothesis. The list accessed on November 2015 had 218 empires with a mean maximum size of 1.161 square miles, with a minimum of 0.060 sq.miles (old kingdom Assyria in 1730 BC) and a maximum of 13 sq.miles (British Empire as measured in the year 1922). We find a significant positive correlation (pearson's correlation coefficient r= 0.29, p-value = 1.329e-05) between the size of empires and the "era" in which they happened reach their zenith ( also see figure below). The correlation is not driven by a few outliers and still shows a strong pattern when restricted to empires smaller than 0.5 million square miles (r=0.16) or 0.1 million square miles (r=0.39).


This suggests that empires have obviously become larger over time. It can also be seen that the largest of them all, the British Empire covered ~23% of the world's land area. One can wonder (when modelling the process is beyond..), if this is some sort of a limit on how big an empire can get before running into trouble. Keith Jeffrey in his book "The British Army and the Crisis of Empire, 1918-22: page 160", echoes this very thought expressed by 'Boney' Fuller when he said "we cannot increase the size of the army to fit the Empire; consequently there is only one thing we can do, namely, reduce the size of the Empire to fit our army".

The age of empires has interesting patterns that are potentially driven by complex processes. Would this correlation exist on a planet that completely lacked geographical barriers? How different would this pattern be if technology had not grown as fast as it has on our blue sphere?


Monday, April 6, 2015

How common is positive selection "inference" in different parts of the human genome?

The database of positive selection provides a handy resource for analyzing patterns of population specific positive selection inference trends and their prevalence in different parts of the human genome. After converting the database into a bed file (loosing some records during conversion makes it easier to run a bunch of analysis on the dataset. 

We next divide the genome into 50Kb chunks and count number of distinct populations in which positive selection was inferred in that window. By this method we can identify "hot-spots" of positive selection in the human genome. Below plot shows the number of populations in which positive selection inference has been done in each of these windows. Only one window on chromosome 2 stretching from 72500000 to 72550000 has been implicated in more than 20 populations. While the Y-chromosome and Mt are empty :) the X chromosome has ~2271 windows with atleast one population.

The above plot shows that it some regions of the genome have been implicated in multiple populations more often than others. Such regions could be termed "hot-spots" of positive selection. However, it reflects the bias in our ability to identify/infer such events rather than actual differences in the occurrence of such events.

Are some populations having positive selection inference spread out over the genome than others? Is it possible to make such a comparison with heterogeneity in the dataset? Methods used are different as well as the resolution of the scans. However, we see(in below figure) how the distribution of such inferences look for populations that have spanned more than 1000 of the ~62,000 windows (50 Kb) spanning the human genome.


MKK has the most number of windows followed by CHB, CEU and YRI. Since, not all observations are independent, it is hard to draw any broadscale conclusions from these numbers.