Mercurial > repos > ucsb-phylogenetics > osiris_phylogenetics
view 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 source
#!/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";