How to simulate NGS reads, controlling sequence coverage?
I am working on a Illumina sequencing simulator for metagenomics: InSilicoSeq

It is still in alpha release and very experimental, but given a multi-fasta and an abundance file, it will generate reads from your input genomes with different coverages.

From the documentation:

iss generate --genomes genomes.fasta --abundance abundance_file.txt \
    --model_file HiSeq2500 --output HiSeq_reads


Where:

# multi-fasta file
&gt;genome_A
ATGC...
&gt;genome_B
CCGT...
...

# abundance file (total abundance must be 1!)
genome_A    0.2
genome_B    0.4
...


I didn't design it to work with coverage but rather abundance of the genome in a metagenome, so you might have to do a tiny bit of math ;)