How To Visualize Splice Junction Data, In Particular Exon-Skipping

Question

[bruce.moran]

Hi all,

I am looking for a method to graph splice junction data, in particular exon-skipping. I have exon-junctions per transcript in count form and just want to be able to graph this to show potential isoforms. Splicegrapher looks nice but is impenetrable to me due to a lack of *nix/python. I have tried GenomicFeatures for plotting the transcripts but cannot find a way to put in an 'alternative splice event' (ie a line connecting exons and missing out those exons for which there is evidence of skipping).

I have R and PERL, so trying to keep to these, or simpler 'packages'. Looking for something like the front page of Splicegrapher.

All ideas welcome!

[Whetting]

I know it doesn't really fix the problem, but how about using cygwin? You should be able to run your program that way. The manual to Splicegrapher looks pretty detailed, so not knowing python should not be a (big) problem.

[bruce.moran]

I use cygwin locally, I'm ok at Unix and use a server for computationally intensive stuff, its just when I get errors out of Splicegrapher or FDM or other python based things I dont have a clue what to do to fix them. My commands seem ok (to me!).

The manual is indeed detailed, but when I copy/paste the commands it gives an error. So I slightly change the command to what I hope might work but still, same error. Eg for FDM (also quite well 'annotated'):

#manual says:

 ./fdm_main.py --config <root>/runs/run2600/fdm_2600.cfg

#when I run this from my FDM_dist dir:

 ./fdm_main.py --config /home/my/stuff/runs/run2600/fdm_2600.cfg
 -bash: root: No such file or directory

#fair enough, obviously meant to change <root> to /home/my/stuff/

 /home/my/stuff/FDM_dist/src/fdm_main.py --config /home/my/stuff/runs/run2600/fdm_2600.cfg
 -bash: ./fdm_main.py: No such file or directory

#ok, lets find where that .py is and make the command point to it

 Traceback (most recent call last):
 File "/home/my/stuff/FDM_dist/src/fdm_main.py", line 193, in <module>
 param_dict=verifyfdmconfig(cfgfile)
 File "/home/my/stuff/FDM_dist/src/fdm_main.py", line 94, in verifyfdmconfig
 root_dir = config.get('Directory','root_dir')
 File "/usr/lib64/python2.6/ConfigParser.py", line 532, in get
 raise NoSectionError(section)
 ConfigParser.NoSectionError: No section: 'Directory'

And then I look at that mess and decide I dont like python and want it in R or PERL please=)

Answer 1

[Julien Textoris]

Hi,

maybe you could use annmap. It allows to plot various kind of numerical data on top of genomic representation. I usd i once with exon microarray data, and plotting the level of expression, according to chromosome coordinate just made it very easy to see which exon was expressed or not.

You will need for this the libraries xmapbridge and xmapcore. I remember the installation was tricky but well documented (you need ensembl data if i remenber well).

Below is a piece of my old script (don't know if it still works, but it may help you to start ?) :

library(xmapbridge)
library(xmapcore)
xmap.connect()

coord.match.best = read.table("/data/tmp/annotation_files/MoExVsMoGene_BestMatch.txt",
                          sep="\t",quote="",header=T)
coord.match.good = read.table("/data/tmp/annotation_files/MoExVsMoGene_GoodMatch.txt",
                          sep="\t",quote="",header=T)
coord.match.complex = read.table("/data/tmp/annotation_files/MoExVsMoGene_Complex.txt",
                          sep="\t",quote="",header=T)
coord.match.no = read.table("/data/tmp/annotation_files/MoExVsMoGene_NoMatch.txt",
                          sep="\t",quote="",header=T)

colnames(coord.match.best)

g = symbol.to.gene("Nfkb1",as.vector=FALSE)

pst = gene.to.probeset(g[["stable_id"]],as.vector=F)

g.probes.MoGe = as.vector(na.omit(names(y[y%in%"Nfkb1"])))

as.character(coord.match.best$A.Probe.Set.Name[which(
        as.character(coord.match.best$B.Probe.Set.Name)%in%g.probes.MoGe)])
as.character(coord.match.good$A.Probe.Set.Name[which(
        as.character(coord.match.good$B.Probe.Set.Name)%in%g.probes.MoGe)])
as.character(coord.match.complex$A.Probe.Set.Name[which(
        as.character(coord.match.complex$B.Probe.Set.Name)%in%g.probes.MoGe)])
as.charactercoord.match.no$A.Probe.Set.Name[which(
        as.charactercoord.match.no$B.Probe.Set.Name)%in%g.probes.MoGe)])

chr = as.vector(unlist(g[["space"]])) 

start = min(as.numeric(as.vector(unlist(g[["ranges"]]))))
end = max(as.numeric(as.vector(unlist(g[["ranges"]]))))

annots.probes=read.table("/data/tmp/annotation_files/MoGene-1_0-st-v1.na32.mm9.probeset.csv",
                     sep=",",quote="\"",header=T)
colnames(annots.probes)

probes = as.vector(na.omit(names(y[y%in%"Nfkb1"])))
starts = annots.probes$start[which(annots.probes$probeset_id%in%probes)]
ends = annots.probes$stop[which(annots.probes$probeset_id%in%probes)]

xc = starts+(ends-starts)/2  ## Here you have the chromosomal coordinate (the x axis)
yc = round(2^fc4nfkb_subunits[which(row.names(fc4nfkb_subunits)%in%probes),1],2)  ## Here you have my expression values (y axis)

xmap.plot(xc,yc,chr,species="mus_musculus",type="area",
      col="red",ylim=c(0,16))

yc2 = round(2^fc4nfkb_subunits[which(row.names(fc4nfkb_subunits)%in%probes),4],2)
xmap.points(xc,yc2,chr,type="area",col="blue")

yc3 = round(2^fc4nfkb_subunits[which(row.names(fc4nfkb_subunits)%in%probes),7],2)
xmap.points(xc,yc3,chr,type="area",col="green")

[bruce.moran]

That looks similar to a function in GenomicFeatures I have used to graph transcripts, GenomeGraphs. It is appreciated, but I need something that can show, on one graph, exon-skipping a la SpliceGrapher or FDM: http://splicegrapher.sourceforge.net/_images/AT5G22640_pred.png