	Computational and Analytical Molecular Evolution Lab at CARB						Nexplorer
	People	Publications	Software	Opportunities	Links	Internal	Nexplorer

NEXUS Stuff: Perl API, demo tools, test files, etc.

Please go to our separate NEXUS projects page to find out more about Bio::NEXUS, our NEXUS applications programming interface.
The Bio::NEXUS package at CPAN includes the library and the nexplot and nextool applications.
If you wish to explore the functionality of nexplot without installing the library, please try the Nexplorer server.
For those of you developing or testing NEXUS software, the NEXUS projects page has information and links for our NEXUS test-file set.

Miscellaneous applications software

fasfil.pl, a handy filter for dealing with fasta archives

BDM.pl, original Perl prototype of the Bateson-Dobzhansky-Muller model used in the Yampolsky-Stoltzfus paper to demonstrate the effect of a bias in the introduction of variation

mb-coda, scripts to evaluate whether the Markov chains of MrBayes have converged

mol2con.pl, a simple script used in making consensus trees from molphy output

ABaCUS, Analysis of Blake's Conjecture Using Simulations

pip, a tool for parsing intron position data from GenBank files

generic.pl, provides command-line interface & help system, very convenient for writing scripts

The development of Intrez, an Entrez client for collating intron data on a gene family, has been dropped.

Stupid Excel Tricks

For didactic reasons, its useful to implement algorithms with software that almost everyone has and almost anyone can use. For better or worse, that means using Excel. Note that you cannot read the files below with your browser-- instead, save the file to disk and open it with Excel (or with StarOffice or gnumeric on your UN*X workstation).

This one-dimensional cellular automaton makes the kind of fractal triangle pattern seen on some snail shells (warning: color plot will take many seconds to re-compute each time you click on it)
This pairwise alignment spreadsheet shows how match scores and path scores are computed in the Needleman-Wunsch (or Smith-Waterman) algorithm. The accompanying explanation of the underlying Excel commands will be of interest to any brave geek who wishes to extend the algorithm to include gap penalties or an arbitrary matrix of match scores (any takers?).

R tutorial

The main R web site is an excellent resource.
A good tutorial exercise is in Appendix A of An Introduction to R

My tutorial exercises are given below (or get the commands in plain text from tutorial.R)

simple examples of an array made with the c (concatenate) command

x1 = c(0.2, 0.4, 0.6, 0.8, 1.0)
x1
 [1] 0.2 0.4 0.6 0.8 1.0
x2 = (1:5)/5
x1 == x2
 [1] TRUE TRUE TRUE TRUE TRUE

here's how we can use an "index vector" to extract the part of a vector we want

index = x1 < 0.5
index
 [1] TRUE TRUE FALSE FALSE FALSE
x1[index] 
 [1] 0.2 0.4

now lets create some x and (randomized) y values, plot them, and add a line with abline

x = (1:100)/5
y = x + rnorm(n=length(x), mean=0, sd=1)
plot(x,y)
abline(lm(y ~ x), col="red")

same thing, but with a sine thrown in to create phasic data. We'll plot smoothed data, then add the data points using the points command.
```
x = (1:100)/5
y = rnorm(length(x),x+sin(x),1)
plot(smooth.spline(x,y), type="l", col="blue")
points(x,y)
```

here we'll load vote data using read.table, and use text to show points with strings (file is fl2000.txt)

florida = read.table("fl2000.txt",header=TRUE)
totals = apply(florida[4:17],1,sum) # compute total ballots per county
florida = cbind(florida,totals)
rm(totals)
plot(florida$totals,florida$Buchanan)
# use type="n" (no points) so we can plot county names in florida$county
plot(florida$totals,florida$Buchanan, type="n", xlab="Total ballots for county", ylab="Votes for Buchanan")
text(florida$totals,florida$Buchanan, florida$county)
abline(lm(florida$Buchanan ~ florida$totals),col="red")

now we'll apply what we've learned to a mass-spec example (input file hyal1iso1f.txt)

read in file

ms.data <- read.table("hyal1iso1f.txt",header=TRUE)
dim(ms.data)  # OOPs, got only one column of data-- I need to specify separator=tab
ms.data <- read.table("hyal1iso1f.txt",header=TRUE, sep="\t")
dim(ms.data) # correct result of 42000+ rows and 2 columns
ms.data[1,]  # see the first row
ms.data[1:10,]  # see the first 10 rows

plot the data

plot(ms.data,pch=".")  # plot all of the data as tiny points

use index vectors based on the logical AND (&) to grab the subset of the data where the peak is. Then we will plot a SMOOTHED version.

peakregion <- ms.data[ms.data[,1]>45000 & ms.data[,1]<57000,]
plot(smooth.spline(peakregion), type="l", col="blue", xlab="m/z", ylab="intensity")
points(peakregion,pch=".")  # good scientists always show the real data, not just the smoothed (fake!) results

make a nice text label for the peak

intensity.max.y = max(peakregion[,2])
intensity.max.x = peakregion[peakregion[,2]==intensity.max.y,1]
peak.label = paste("max at ", format(intensity.max.x,nsmall=1))
text(intensity.max.x+3000, intensity.max.y-200, peak.label)

finally, we use plot, points, text and get output in PostScript

postscript(file="ms.data.plot.ps",horizontal=FALSE, width=5,height=5)
plot(smooth.spline(peakregion,nknots=50),type="l", col="blue",xlab="m/z",ylab="intensity",main="Main peak region")
points(peakregion,pch=".")
text(intensity.max.x+3000, intensity.max.y-200, peak.label)
dev.off()  
# "dev.off()" is needed because for R, PostScript is a live stream to an output
# device-- we have to shut off the tap now that we are done with it