Mercurial > repos > ucsb-phylogenetics > osiris_phylogenetics
diff phylographics/makeRtrees.pl @ 0:5b9a38ec4a39 draft default tip
First commit of old repositories
author | osiris_phylogenetics <ucsb_phylogenetics@lifesci.ucsb.edu> |
---|---|
date | Tue, 11 Mar 2014 12:19:13 -0700 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/phylographics/makeRtrees.pl Tue Mar 11 12:19:13 2014 -0700 @@ -0,0 +1,121 @@ +#!/usr/bin/perl + +#This script generates an R script to print trees to a pdf file +#input is a table with treename<tab>newick tree +use strict; + +my $filename = $ARGV[0]; +my $outfile = $ARGV[1]; +open FILE, $filename or die $!; +my $treetype = $ARGV[2]; +my $extiplabels = $ARGV[3]; +my $options; +my $labeltaxfile = $ARGV[4]; +my $query = $ARGV[5]; +my $midroot = $ARGV[6]; +my %labelhash; +my $genecount=0; +my @genes; +my $markquery; +my $midpoint; + +if($query eq 'yes'){ + $markquery = 1; +}else{ + $markquery = 0; +} +unless($labeltaxfile eq 'None'){ + open LABELFILE, $labeltaxfile or die $!; + while (<LABELFILE>) { + chomp; + #get a line from the data file + my $currentinput = "$_"; + if($currentinput =~ /\t/){ + my @splitline = split(/\t/); + my $speciesname= $splitline[0]; + $speciesname = "'".$speciesname."'"; + my $treename = $splitline[1]; + if(exists $labelhash{$treename}){ + push @{ $labelhash{$treename} }, $speciesname; + }else{ + push @{ $labelhash{$treename} }, $speciesname; + #$labelhash{$treename} = $speciesname; + $genecount ++; + push @genes, $treename; + } + } + } + +}#end unless + + +if($extiplabels eq 'yes'){ + $options = ", show.tip.label=FALSE"; +}else{ + $options = ", show.tip.label=TRUE"; +} + +print "require(ape);\n"; +print "require(phangorn);\n"; +print "pdf(file='$outfile');\n"; + +while (<FILE>) { + chomp; + #get a line from the data file + my $currentinput = "$_"; + my @splitline = split(/\t/); + my $treename= $splitline[0]; + my $tree = $splitline[1]; + my $labelsvector; + + #print the R commands to make tree graphics + print "raw_tree <- read.tree(text = '$tree');\n"; + print "raw_tree\$edge.length[ is.na(raw_tree\$edge.length) ] <- 0 \n"; + if($midroot eq 'yes'){ + print "raw_tree <- midpoint(raw_tree)\n"; + } + #Check if large tree, then make text size smaller + print "thetips <- raw_tree\$tip.label \n"; + print "numtips <- length(thetips) \n"; + +#Make text smaller for trees with many tips + print "if(numtips>250){plot(raw_tree, cex=0.15, edge.width = 0.1, type='$treetype' $options)}else if(numtips>75){plot(raw_tree, cex=0.3, type='$treetype' $options)}else{plot(raw_tree, cex=0.6, type='$treetype' $options)};\n"; + print "title('Tree File: $treename');\n"; + +#Add taxon labels, if optional file present and if labels exist for tree + if(exists $labelhash{$treename}){ + $labelsvector = join ",", @{ $labelhash{$treename} }; + $labelsvector = "tolabel <- c(".$labelsvector.")"; + print "thetips <- raw_tree\$tip.label \n"; + print $labelsvector."\n"; + print "labels <- match(tolabel,thetips) \n"; + print "tiplabels(tip=labels, pch=21, cex=1) \n"; + } + + +#Add taxon labels if gene name contains QUERY - for readplacement + if($markquery == 1){ + print "thetips <- raw_tree\$tip.label \n"; + print "qlabels <- grep(\'QUERY\',thetips) \n"; + print "tiplabels(tip=qlabels, pch=21, cex=1.1) \n"; + print "l1labels <- grep(\'LANDMARK1\',thetips) \n"; + print "tiplabels(tip=l1labels, pch=15, cex=.8, col='red') \n"; + } +} + +print "dev.off();\n"; +close FILE; + +#Testing hash arrays +#my %nums; +#my $test='odd'; +#for my $n (4,5,6,10) { +# if ($n % 2) { +# push @{ $nums{$test} }, $n; +# } else { +# push @{ $nums{even} }, $n; +# } +#} +# +#print join ', ', @{ $nums{even} }; +#print "\n\n";