Hi, I would like to ask if you can recommend me tools that I can use to simulate whole genome sequences using a reference genome which will also give me list of the variants just like wgsim?
I was able to use wgsim, but when I tried aligning the paired-reads, it's not aligning together.
reads=100000000
for i in $(seq 1 10); do
base="sim_${i}"
seed=$((100 + i)) # different seed for each run (arbitrary choice)
echo "Running $base (seed=$seed)..."
"$wgsim" "$ref" "${base}.R1.fq" "${base}.R2.fq" -1 70 -2 70 -N $reads -S $seed -e 0.0001 > "${base}.out.log" 2>&1
echo "$base finished (log: ${base}.out.log)"
done
of maybe I am using wgsim wrong? I hope you can help me. Thank you!
3 answers
The best one I've used for ONT reads so far is badread - https://github.com/rrwick/Badread
For illumina I've used and like insilicoseq - https://insilicoseq.readthedocs.io/en/latest/
Wow, thats too long. How much RAM does your machine have, it is likely swapping
You may want to consider one of these packages:
I think only a fraction of sequencing errors can be specified. If a mix of a reference and a mutated genome is included, setting sequencing errors to 0 might give reads with fixed mutations.
if you can recommend me tools
You can use randomreads from BBMap suite to generate short illumina reads with known mutations. A guide is available here: https://bbmap.org/tools/randomreads
Log in to answer this question.
-
mRNA Seq
written by Yaseen •
-
How to run MLST with multiple fastq files
written by Kumar
-
Should I merge the Paired-end RNA-seq reads for BLAT for intron identification
written by praasu •
-
run multiple files for sortmeRNA
written by tianshenbio
-
Obtain Number of PE Reads in BAM file Mapped to Each Chromosome
written by ssm87 •
-
bbduk: quality filtering and singletons
written by rororo •
-
Remove duplicates reads Ids
written by BioGeek
-
understanding wgsim flags
written by kshitijtayal •
-
How to run Spades For Nextseq data
written by jeccy.J
-
How To Detect/Filter Sequences That Contain A Large Number Of N Symbols?
written by PoGibas
It's a little unclear what you mean by "I was able to use wgsim, but when I tried aligning the paired-reads, it's not aligning together."
Sorry about that. It's when I align the paired reads together like the left and righ, they should align with common sequence at some ends, making one longer consensus sequence. I was able to align the pairs with an actual short read pairs.
I believe the pairs will only have common sequence (e.g. the pairs will "overlap") if the insert size is small. the wgsim program has a flag called -d ("outer distance between the two ends") which i think can adjust the insert size and might be able to be made smaller (default: 500) to make it produce some overlap but I haven't tested it myself
You don't want to simulate reads like this. Good WGS libraries should not have reads that overlap in the middle because they will represent short inserts. Not what one wants in real life.
You want more data to cover a particular mutation intriduced by simulation by having more read pairs covering it, rather than having anoverlapping read pair covering that mutation.
I did a quick search for tools in the NIAID Data Ecosystem and found a few that map Illumina short reads onto a reference genome, including Stampy and Novocraft. But you may be able to find something better or more specific using the tool filters.
This question is about simulating reads and not read mapping.