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";