SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Scripture - transcriptome reconstruction Boel Bioinformatics 17 06-27-2012 11:55 AM
RNA-Seq: Full-length transcriptome assembly from RNA-Seq data without a reference gen Newsbot! Literature Watch 7 10-26-2011 06:37 AM
RNA-seq: ab initio reconstruction of transcriptome with prokaryotes pasta RNA Sequencing 3 07-26-2011 08:22 PM
RNA-Seq: Composite Transcriptome Assembly of RNA-seq data in a Sheep Model for Delaye Newsbot! Literature Watch 0 03-26-2011 03:02 AM
RNA-Seq: Accurate quantification of transcriptome from RNA-Seq data by effective leng Newsbot! Literature Watch 0 11-10-2010 03:00 AM

Reply
 
Thread Tools
Old 08-06-2011, 07:19 AM   #1
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Arrow Transcriptome reconstruction from RNA-seq data using Scripture.

Hi all,

I'm trying to re-construct the transcriptome from RNA-seq data using scripture and I get some errors. I'm pasting a log of errors and commands I've been using.

I have paired end data from Illumina, I use tophat on these files as:
Code:
$nohup tophat -o tophat_out_SRR065506_1 --GTF ../ucsc.gtf ../hg19index/hg19-25chr SRR065506_1.fastq

$nohup tophat -o tophat_out_SRR065506_2 --GTF ../ucsc.gtf ../hg19index/hg19-25chr SRR065506_2.fastq
I get the standard output files. The accepted_hits.bam is converted into accepted_hits.sam as:

Code:
$ samtools view -h -o tophat_out_SRR065506_1/accepted_hits.sam tophat_out_SRR065506_1/accepted_hits.bam

$ samtools view -h -o tophat_out_SRR065506_2/accepted_hits.sam tophat_out_SRR065506_2/accepted_hits.bam
I now take off the header as:

Code:
$ sed '1,2d' tophat_out_SRR065506_1/accepted_hits.sam | sort > tophat_out_SRR065506_1/accepted_hits.sorted.sam

$ sed '1,2d' tophat_out_SRR065506_2/accepted_hits.sam | sort > tophat_out_SRR065506_2/accepted_hits.sorted.sam
There is only one apparent change in file. I attach screenshots of a small part of the two files.





I then use Scripture's -task makePairedFile as:

Code:
$java -jar scripture.jar -task makePairedFile -pair1 tophat_out_SRR065506_1/accepted_hits.sorted.sam -pair2 tophat_out_SRR065506_2/accepted_hits.sorted.sam -out postTophat/SRR065506.scripturePaired.sam -sorted
On completion, I use IGVtools to sort and index the paired and alignment files as:
Quote:
$cat tophat_out_SRR065506_1/accepted_hits.sorted.sam tophat_out_SRR065506_2/accepted_hits.sorted.sam > postTophat/all_tophat_alignments.sam
$igvtools sort postTophat/all_tophat_alignments.sam all_tophat_alignments.sorted.sam

$igvtools sort postTophat/SRR065506.scripturePaired.sorted.sam
This completes successfully and I get the .sai files.

I then use Scripture as:

Code:
$ java -jar scripture.jar -alignment all_tophat_alignments.sorted.sam -out scriptureResults/chr1.segment -sizeFile ../hg19/hg19.chrom.sizes2 -chr chr1 -chrSequence ../hg19/chr1.fa -pairedEnd SRR065506.scripturePaired.sorted.sam
This is when I get the error. Log of the error...

Code:
Using Version VPaperR3
Computing weights..... upweighting? false weight: 1.0
Computing alignment global stats for chromosome chr1
Computing alignment global stats for chromosome chr10
Computing alignment global stats for chromosome chr11
Computing alignment global stats for chromosome chr12
Computing alignment global stats for chromosome chr13
Computing alignment global stats for chromosome chr14
Computing alignment global stats for chromosome chr15
Computing alignment global stats for chromosome chr16
Computing alignment global stats for chromosome chr17
Computing alignment global stats for chromosome chr18
Computing alignment global stats for chromosome chr19
Computing alignment global stats for chromosome chr2
Computing alignment global stats for chromosome chr20
Computing alignment global stats for chromosome chr21
Computing alignment global stats for chromosome chr22
Computing alignment global stats for chromosome chr3
Computing alignment global stats for chromosome chr4
Computing alignment global stats for chromosome chr5
Computing alignment global stats for chromosome chr6
Computing alignment global stats for chromosome chr7
Computing alignment global stats for chromosome chr8
Computing alignment global stats for chromosome chr9
Computing alignment global stats for chromosome chrM
Computing alignment global stats for chromosome chrX
Computing alignment global stats for chromosome chrY
Has pairs: true
Has upweighting turned on: false
Computing weights..... upweighting? false weight: 1.0
AlignmentDataModel loaded, initializing model stats
Computing alignment global stats for chromosome chr1
Computing alignment global stats for chromosome chr10
Computing alignment global stats for chromosome chr11
Computing alignment global stats for chromosome chr12
Computing alignment global stats for chromosome chr13
Computing alignment global stats for chromosome chr14
Computing alignment global stats for chromosome chr15
Computing alignment global stats for chromosome chr16
Computing alignment global stats for chromosome chr17
Computing alignment global stats for chromosome chr18
Computing alignment global stats for chromosome chr19
Computing alignment global stats for chromosome chr2
Computing alignment global stats for chromosome chr20
Computing alignment global stats for chromosome chr21
Computing alignment global stats for chromosome chr22

Exception in thread "main" net.sf.samtools.SAMFormatException: Error parsing text SAM file. Not enough fields; Line 510022
Line: @SQ	SN:chr15	LN:102531392
	at net.sf.samtools.SAMTextReader.reportFatalErrorParsingLine(SAMTextReader.java:169)
	at net.sf.samtools.SAMTextReader.access$400(SAMTextReader.java:40)
	at net.sf.samtools.SAMTextReader$RecordIterator.parseLine(SAMTextReader.java:268)
	at net.sf.samtools.SAMTextReader$RecordIterator.next(SAMTextReader.java:232)
	at net.sf.samtools.SAMTextReader$RecordIterator.next(SAMTextReader.java:196)
	at org.broad.igv.sam.reader.SamQueryTextReader$SAMQueryIterator.next(SamQueryTextReader.java:197)
	at org.broad.igv.sam.reader.SamQueryTextReader$SAMQueryIterator.next(SamQueryTextReader.java:141)
	at broad.pda.seq.segmentation.GenericAlignmentDataModel.getCountsPerAlignment(GenericAlignmentDataModel.java:257)
	at broad.pda.seq.segmentation.GenericAlignmentDataModel.getCountsPerAlignment(GenericAlignmentDataModel.java:196)
	at broad.pda.seq.segmentation.GenericAlignmentDataModel.getCountsPerAlignment(GenericAlignmentDataModel.java:1661)
	at broad.pda.seq.segmentation.AlignmentDataModelStats.getNumberOfReadsByChr(AlignmentDataModelStats.java:178)
	at broad.pda.seq.segmentation.AlignmentDataModelStats.computeDataStats(AlignmentDataModelStats.java:133)
	at broad.pda.seq.segmentation.AlignmentDataModelStats.computeGlobalStats(AlignmentDataModelStats.java:122)
	at broad.pda.seq.segmentation.AlignmentDataModelStats.<init>(AlignmentDataModelStats.java:86)
	at broad.pda.seq.segmentation.AlignmentDataModelStats.<init>(AlignmentDataModelStats.java:67)
	at broad.pda.seq.segmentation.ContinuousDataAlignmentModel.main(ContinuousDataAlignmentModel.java:2277)
I don't reckon that I've made a major mistake while doing this. My most probable guess is that something is wrong in format conversion (bam to sam). I'm working on that right now.

Any insights or help would be greatly appreciated.

Thanks!

Last edited by Joker!sAce; 08-06-2011 at 07:22 AM.
Joker!sAce is offline   Reply With Quote
Old 08-06-2011, 09:43 AM   #2
upendra_35
Senior Member
 
Location: USA

Join Date: Apr 2010
Posts: 102
Default

Hi

I have exactly got the same error as you "Exception in thread "main" net.sf.samtools.SAMFormatException: Error parsing text SAM file. Not enough fields". I don't know where i am doing wrong. Let me know if you find the solution for this problem. However there is another way of doing scriputre analysis which i did and it worked. Just follow the scripture manual in this link (http://www.molecularevolution.org/re...pture_activity). Good luck!
upendra_35 is offline   Reply With Quote
Old 08-07-2011, 03:26 AM   #3
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

Hi,

Thanks for the reply.

I have 2 questions on this link...

1. In the IGVtools count command, what is the .genome file? Rather how do you get it? I tried the .fa reference genome but that did not work.
2. In scripture, makePairedFile task requires 2 files. In the link, thye've done it in a different way. I unzipped the scripture and scripture_beta versions but I cannot find the file scripture_alpha2.jar What/Where is this file?

Thanks!

Last edited by Joker!sAce; 08-07-2011 at 03:33 AM.
Joker!sAce is offline   Reply With Quote
Old 08-07-2011, 03:38 AM   #4
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

I was using IGVtools without the genome files. Downloading that, it should hopefully work now. Q2 still stands.

Last edited by Joker!sAce; 08-07-2011 at 03:52 AM.
Joker!sAce is offline   Reply With Quote
Old 08-07-2011, 01:33 PM   #5
upendra_35
Senior Member
 
Location: USA

Join Date: Apr 2010
Posts: 102
Default

Hi,
I guess the workshop material for Scripture analysis was run on old version (alpha) and the current downloaded version for broad institute was new (beta). I think you can use the beta version and see if it still works or otherwise i will send you the old version.
upendra_35 is offline   Reply With Quote
Old 08-07-2011, 02:49 PM   #6
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

Hey,

That'd be awesome. My e-mail is orsha.pk@gmail.com

I've tried the three versions of scripture they've put on the website, none of them have this older function. Maybe they should include the older functions too in the next release, I'm going to leave a comment there.

Thanks!
Joker!sAce is offline   Reply With Quote
Old 08-07-2011, 06:41 PM   #7
chenyao
Member
 
Location: Beijing

Join Date: Jul 2011
Posts: 74
Default

why don't you use cufflinks?
chenyao is offline   Reply With Quote
Old 08-08-2011, 02:30 AM   #8
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

Everyone in my department already does. We wanted to test out scripture. I know that using cufflinks would probably be less challenging given the great documentation.

Last edited by Joker!sAce; 08-08-2011 at 02:32 AM.
Joker!sAce is offline   Reply With Quote
Old 08-10-2011, 02:10 AM   #9
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

@ Upendra

Thanks a lot for the mail. It all works well now.

I have one small question.

In the bed file there are some transcripts with a * instead of +/- for orientation. What does this mean?
Joker!sAce is offline   Reply With Quote
Old 08-16-2011, 06:58 AM   #10
chenyao
Member
 
Location: Beijing

Join Date: Jul 2011
Posts: 74
Default

what's the difference between scripture and cufflinks,

I also have a question, scripture reconstruct transcriptome by specific chromosome, How can I use whole genome? I don't see the command. And where is the whole genome size file?
chenyao is offline   Reply With Quote
Old 08-17-2011, 05:44 AM   #11
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

A size file is just a file with chromosome name and it's length. i've attached 2 text files. The chrom2 is the one you should probably use, the other has lengths of all hg19 chromosomes(alternative haplotypes).

Difference in scripture and cufflinks is in how it reconstructs the transcriptomes. Read the papers linked on the official sites. Also, there is plus and minus of each software. While cufflinks can handle whole genome, differential analysis, basic stats; it has aggressive filtering (and more that I don't know of). Scripture does not have that aggressive filtering, but scripture pipelines are a bit difficult to manage.

You cannot run scripture for the whole genome. What you can do it make a bash script that runs for all the chromosomes, that'd automate running the same command over and over again. Concatenate all the .bed files. Sort it and covert to .gtf file.
Attached Files
File Type: txt hg19_chrom_sizes2.txt (375 Bytes, 21 views)
File Type: txt hg19_chrom_sizes.txt (1.9 KB, 9 views)

Last edited by Joker!sAce; 08-17-2011 at 05:47 AM.
Joker!sAce is offline   Reply With Quote
Old 01-28-2012, 06:03 AM   #12
nbartonicek
Junior Member
 
Location: Cambridge, UK

Join Date: Jul 2011
Posts: 3
Default

The issue from the first post was caused by the creation of the sam files with the "-h" option. This did not allow the proper cleaning of the headers, which resulted in erroneous merged files.
...
Talking from experience here.

Cheers,

Nenad
nbartonicek is offline   Reply With Quote
Old 07-04-2012, 10:33 AM   #13
oxydeepu
Member
 
Location: bangalore,india

Join Date: Jul 2011
Posts: 41
Default

Hi all,

Am following the guide lines mentioned in the scripture walk through example.
But the command

java -Xmx4000m -jar scripture.jar -task makePairedFile -pair1 tophat_out_SRR039999_1/accepted_hits.sorted.sam -pair2 tophat_out_SRR039999_2/accepted_hits.sorted.sam -out SRR039999.paired.sam -sorted

This is not giving me any kind of out put can anyone help me out.

Thank you in advance
Deepak
oxydeepu is offline   Reply With Quote
Old 12-06-2013, 12:31 PM   #14
11xinqi
Member
 
Location: china

Join Date: Mar 2011
Posts: 31
Default

Quote:
Originally Posted by upendra_35 View Post
Hi,
I guess the workshop material for Scripture analysis was run on old version (alpha) and the current downloaded version for broad institute was new (beta). I think you can use the beta version and see if it still works or otherwise i will send you the old version.
Hi, I want to use the scripture (old version) to transfer bed file to gff file. BC the beta version can not do that. My email is joseph_wangj@hotmail.com Would you please send me a copy of old version
11xinqi is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 01:12 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO