RNAseq for DE purpose

Hi all,

I am totally new in the bioinformatic analysis. I am working on a project that looks at DGE among different time treatments. Besides, there is no reference genome (meaning that I need a de novo assembly step).

So far, after struggling and navigating many website, I started using my terminal on MacBook and was able to do FastQC for the 118 reads (PE), I then followed this with MultiQC and both worked fine for me using the terminal. I am trying to go one step at a time. So I wanted to use trimmomatic to trim the reads but I failed!

I used the following command for one pair of reads and it worked but I do not know how to do this for the whole reads!

java -jar trimmomatic-0.39.jar PE -phred33 $i\1_R1_001.fastq.gz $i\1_R2_001.fastq.gz $i\1_R1_paired.fq.gz $i\1_R1_unpaired.fq.gz $i\1_R2_paired.fq.gz $i\1_R2_unpaired.fq.gz ILLUMINACLIP:contams_forward_rev.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36

I tried adding * but it did not work.

java -jar trimmomatic-0.39.jar PE -phred33 -threads 4 /Users/mohamedhabib/Desktop/RNAseq/*_R1_001.fastq.gz *_R1_001.fastq.gz out_*_R1_001_P.fastq.gz out_*_R1_001_U.fastq.gz out_*_R2_001_P.fastq.gz *_R2_001_U.fastq.gz ILLUMINACLIP:TruSeq3-PE.fa:2:40:15

I will appreciate your help and advice on this specific topic and generally for further steps of analysis. Is there a template script (containing all the steps) that I can modify and use for my analysis?

Thanks

rnaseq

• 2,941 views • link

updated 3.5 years ago by swbarnes2

• written 3.5 years ago by m.habib •

There is something weird about your commands: what are the \s doing after the variable names?

java -jar trimmomatic-0.39.jar PE -phred33 $i\1_R1_001.fastq.gz $i\1_...
#                                            ^                    ^

• 0 views • link

written 3.5 years ago by Ram

I just copied it as it is from another thread and modified it to my own data

• 0 views • link

written 3.5 years ago by m.habib •

What thread did you copy it from, what changes did you make and what was the logic behind each of those changes?

• 0 views • link

written 3.5 years ago by Ram

here: Trimmomatic job script to run on multiple pair end read file

• 0 views • link

written 3.5 years ago by m.habib •

That answers 1 out of the 3 questions I asked.

• 0 views • link

written 3.5 years ago by Ram

I basically changed the file names to make sure this command is working. I picked one forward and reverse read pair and processed it as below: it worked but when I tried to trim all the samples (I thought maybe by adding a * to catch all the files!!) but it failed.

java -jar trimmomatic-0.39.jar PE -phred33 $i\1_R1_001.fastq.gz $i\1_R2_001.fastq.gz $i\1_R1_paired.fq.gz $i\1_R1_unpaired.fq.gz $i\1_R2_paired.fq.gz $i\1_R2_unpaired.fq.gz ILLUMINACLIP:contams_forward_rev.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36

I hope that make sense because I am not an expert!

• 0 views • link

written 3.5 years ago by m.habib •

You need your command to run once per pair of files, not combine all pairs and then run just once. Read the thread you're trying to emulate and understand the loop part of it. Modify your code after you understand the loop.

• 0 views • link

written 3.5 years ago by Ram

I assume you are referring to 118 samlpes (and not reads). You can use a for loop to do this. See example here: Trimmomatic job script to run on multiple pair end read file

• 0 views • link

written 3.5 years ago by GenoMax

I used the following:

java -jar trimmomatic-0.39.jar PE -phred33 -threads 4 $i\_R1_001.fastq.gz $i\_R2_001.fastq.gz $i\_R1_001¬_paired.fastq.gz $i\_R1_001_unpaired.fastq.gz $i\_R2_001_paired.fastq.gz $i\_R2_001_unpaired.fastq.gz ILLUMINACLIP: TruSeq3-PE.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36 >> cmd_file

and I got:

Exception in thread "main" java.lang.ArrayIndexOutOfBoundsException: Index 1 out of bounds for length 1

• 0 views • link

updated 3.5 years ago by Ram

• written 3.5 years ago by m.habib •

What names do your fastq files have? That first part of my answer will need to be adjusted to account for the differences in your file names so the for loop works right.

• 0 views • link

written 3.5 years ago by GenoMax

Here is an example: 19144R-03-01_S75_L004_R2_001.fastq.gz

• 0 views • link

written 3.5 years ago by m.habib •

Use following:

for i in `ls -1 *R1*.fastq.gz | sed 's/\_R1.fastq.gz//'` do echo trimmomatic PE -phred33 $i\_R1.fastq.gz $i\_R2.fastq.gz $i\_R1_paired.fq.gz $i\_R1_unpaired.fq.gz $i\_R2_paired.fq.gz $i\_R2_unpaired.fq.gz ILLUMINACLIP:contams_forward_rev.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36 >> cmd_file; done

If the commands look correct then remove echo and >> cmd_file parts from the command to actually run the trimming.

• 0 views • link

written 3.5 years ago by GenoMax

Honestly, I think most RNASeq aligners are fairly robust to adapters. You could probably get away with not trimming at all

• 0 views • link

written 3.5 years ago by swbarnes2

0 answers

No answers yet.