Hello all,

I need to sample reads from a concatenation of two genomes to simulate contamination. I have been trying to use seqtk for this, but results have been very confusing. I have a bash script for the batch system of the server, which calls my python script, which in turn calls seqtk for each input FASTQ. Every time I run the script with the same input file set, the output is different. By different, I mean that the number of reads are varying even though it is supposed to be one million across the board. Sometimes all output files are blank. Furthermore, in most runs the number drops as you go down the list of files, e.g.:

4000000 File1

4000000 File2

3943234 File3

3523402 File4

......

50234 File37

24323 File38

0 File39

0 File40

This pattern over multiple runs of seqtk is the most peculiar for me. In my previous uses of seqtk, just by trying over and over I was able to get all complete files. I thought it was a memory issue. This time I am running on the supercluster with 256 GB RAM allocated. I have a bash script for the batch system of the server, which calls my python script, which in turn calls seqtk for each input FASTQ. The batch run completes with a few illegal character errors for 'M', 'R' and 'S'.

The bash script is essentially:

module load biokit

python /path/to/myscript/

used_slurm_resources.bash

The python script is:

#!/usr/bin/env python

import subprocess
import os

seqs = ['list', 'of', 'fastq', 'names']

seedlist = ['3', '16', '43', '50', '59', '69', '73', '74', '100', '102']

args = ['/path/to/seqtk', 'sample', '-2', '-s', '', '1000000']

for i in seqs:
    for j in range(0, 10):
        args[3] = '-s%s' % seedlist[j]
        args[4] = '/path/to/input/%s' % i
        output = '/path/to/output/%ss%d%s' % (i[:20], j, i[19:])
        sample = open(output, 'w')
        seqtk = subprocess.Popen(args, bufsize=-1, stdout=sample)
        sample.close()

Errors I get:

(standard_in) 1: illegal character: M
(standard_in) 1: illegal character: R
(standard_in) 1: illegal character: S
(standard_in) 1: illegal character: S
(standard_in) 2: syntax error
(standard_in) 1: syntax error

EDIT: seqtk seems to be working when I run directly for a single file.