1
$\begingroup$

I have tried to use flux simulation tool to generate simulated RNA-seq data.

I gave the following parameter file to flux-simulation shell script

## File locations

REF_FILE_NAME   cElegansAnnotation.gtf
GEN_DIR         chromFa

## Library preparation
# Expression
NB_MOLECULES    5000000
TSS_MEAN    50
POLYA_SCALE     100
POLYA_SHAPE     1.5

# Fragmentation
FRAG_SUBSTRATE  RNA
FRAG_METHOD UR
FRAG_UR_ETA     350

# Reverse Transcription
RTRANSCRIPTION  YES
RT_MOTIF    default

# Amplification 
PCR_DISTRIBUTION default
GC_MEAN      NaN
PCR_PROBABILITY  0.05

# Size Filtering
FILTERING   NO


## Sequencing
READ_NUMBER 1000000
READ_LENGTH 100
PAIRED_END  YES

# create a fastq file
FASTA           YES

According to this parameter, flux-simulation should have given reads with length 100. However when I look at the output, I have seen that reads with less then length of 100 are also included in my fasta file:

>chrI:47472-49416W:Y48G1C.12:3:651:351:594:A/1
CGTCGAAATTAGTGATATTTTTATCGGGAATCGGTCCGTGTGGTTCTCCGGTGAATATTCGATTCGTTGTGGAGACACGAGATCGCTGGGGTCCAAGGAC
>chrI:47472-49416W:Y48G1C.12:4:651:394:467:S/1
TACGCGACAAAAATGGGAAACCGAATCGCGTTTTTTGGCTTCAAGTACAAGTTATTCAGAATCATCAAAATGGG
>chrI:47472-49416W:Y48G1C.12:4:651:394:467:A/2
CCCATTTTGATGATTCTGAATAACTTGTACTTGAAGCCAAAAAACGCGATTCGGTTTCCCATTTTTGTCGCGTA
>chrI:47472-49416W:Y48G1C.12:5:651:126:142:S/2
AGTTGTAAAAGCGGATT

So, is there anyone how could I fix it?

$\endgroup$
  • $\begingroup$ I think it gives a distribution of read lengths. Not all reads would be 100nt. Even in actual experiments you observe that. Can you add a graph of read length distribution. You can use histogram in matlab, python or R. $\endgroup$ – WYSIWYG May 21 '15 at 5:23
  • $\begingroup$ @WYSIWYG thanks for reply, I'll add it tonight, actually, I could handle with different size of reads but they have to be longer than a specific length. Could it be possible? $\endgroup$ – stackunderflow May 22 '15 at 11:30
  • $\begingroup$ You can remove the small reads before mapping. fastx_toolkit has that option or otherwise just do this: awk 'NR%2==1{h=$0} (NR%2==2) && (length($0)>x){print h"\n"$0}' fastafile.fa $\endgroup$ – WYSIWYG May 22 '15 at 13:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.