Hi all, I’m trying to run Kraken2 in a parallelize...
# help
c
Hi all, I’m trying to run Kraken2 in a parallelized way under a single SLURM job using shared memory mapping (so multiple samples can run concurrently without duplicating large RAM usage). Ideally, I want to process ~128 samples within one SLURM job submitted by nextflow. so far, each input is a tuple:
(sample_ID, sample_fq)
and the expected output is:
(sample_ID, unclassified_fq)
My main issue is how to correctly structure the workflow so that: • input list are automatically split into around 128 samples per list • each list of samples are processed in parallel within the same task/job • each sample’s outputs are correctly tracked and associated with its sample_ID • I avoid launching one SLURM job per sample but still maintain clean sample-level outputs What is the recommended way to design this (e.g. in Nextflow / GNU parallel / custom scripting), and how should I manage input/output bookkeeping for each sample in a single parallelized job? I tried to input and output converting from channel to a list (.csv) in inputing and outputing but I am afraid that this will not resume properly. I have been multithreading to the maximum (128) since it requires the whole node to run large-mem kraken2 job. Or is there any other way apart from parallelizing to save core hours usage on cluster?
n
Your first question could probably be solved with
buffer()
. I think the other things you are trying to do are more complex, because they go a bit against nextflows one sample - one job approach. You would probably need to find a way via parallel or xargs or something that would allow you to then generate the command to launch kraken2 correctly / several times in parallel (not familiar with kraken2)? I would probably then collect all outputs into a single emit (probably named by
sample_id
) and have some channel stuff (combination of
flatten()
and
map()
maybe) to turn them back into a proper channel. The main issue (to me) is probably to figure out if the shared memory thing in slurm actually works as intended in the nextflow generated sbatch script.
❤️ 1
j
We've had the same challenge with KRakenUniq in #C031QH57DSS
c
I have a script to use parallel for many samples in memory mapping in kraken2, tested and works so this is solved. The issue is indeed the attempt to process and keep track of many sample in one job. The number of samples are indefinite so I might need to figure out a way to read and output in case I manage the input as [ [sample_ID, sample], [sample_ID, sample], [sample_ID, sample] ] and output as [ [sample_ID, out_sample], [sample_ID, out_sample], [sample_ID, out_sample] ] and report error if output is missing for a sample_ID
j
(looking for the code, one moment)
❤️ 1
And requires a non-standard module structure
:
But basically we found that you could not really have a shared database ACROSS jobs
We ended up having a single job and loop through it 😕
Which is suboptimal
but the only way around it
c
Thank you James. As long as it works 👍🏻 Is meta a groovy dictionary/map? I saw that it's input as val() but called in script. Could you point me also how it's defined?
Copy code
classified = meta.single_end ? "\${PREFIX}.classified.${sequence_type}" : "\${PREFIX}.merged.classified.${sequence_type}"
    unclassified = meta.single_end ? "\${PREFIX}.unclassified.${sequence_type}" : "\${PREFIX}.merged.unclassified.${sequence_type}"
I am also interested in seeing how the output val(meta), path(unclassified) etc., could be parsed to val(sample_id), path(unclassified), though I did not find them in the shared main.nf. Could you point me where this workflow is called and its output parsed?
j
Yes, correct - it's a standard nf-core style groovy map: https://nf-co.re/docs/developing/components/meta-map
I am also interested in seeing how the output val(meta), path(unclassified) etc., could be parsed to val(sample_id), path(unclassified), though I did not find them in the shared main.nf. Could you point me where this workflow is called and its output parsed?
I'm not sure what you mean exactly here, but are you talking about within the module?
c
I think now I understand. The output does not need to be parsed (so my question dissolves), as previously I thought that meta contain the association between sample_id and sample fq.
j
It does iniitally yes
But then we combine it all into one
But what @Niklas Schandry remains valid
n
Since this may not be obvious: You can stuff pretty much whatever into the meta-map, it is wildly undefined (except that it has to contain
id
). This can be kind of useful if you need to pass things through a process, i.e. you need something afterwards that does not match the process input, you can simply store it somewhere in meta and map it back out later. This can also avoid a lot of joins etc
this 1
❤️ 1
c
Gotcha. I will try with the batching and then meta-map and input list
Thank you two!!
I managed to have many samples in the same process run by using a list
.collectFile
from the previous process output, then batching, and run through kraken2 to output a list which contains the sample_ID and output_files in each row. However, managing the samples with a list than nf channel, now it does not resume KRAKEN2 any more despite a successful run. Since the kraken2 process is at the very beginning of my workflow, not being able to resume from here will break the resumability of the majority of the workflow. Is there any way to restore the resume? Or is it not possible? here are the codes:
Copy code
preprocessed_reads
				.map { it.join(',') }
				.collectFile(name: "${workflow.runName}.csv", newLine: true)
				.set { input_kraken2 }

				input_kraken2.view { "input_kraken2: $it" }

				//split the input of kraken2 into 128 samples per batch
				input_kraken2
				.flatMap { csv_file ->
					// Read the file content and process it
					def lines = csv_file.text.split('\n').findAll { it.trim() } // Remove empty lines

					// Create batches of size 2 (batches of lines, not parsed CSV)
					def batches = lines.collate(5)

					// Return one item per batch
					return batches.withIndex().collect { batch, index ->
						def basename = csv_file.getBaseName()
						def filename = "${basename}.sub${index + 1}.csv"
						def content = batch.join('\n') + '\n'

						return [filename, content]
					}
				}
				.collectFile() { filename, content ->
					[filename, content]
				}
				.set{input_kraken2_subsets}


				input_kraken2_subsets.view { "1_subsets: $it" }

				input_kraken2_subsets
				.map{it}
				.map { it -> tuple(it, params.KRAKEN2_FILTER_DATABASE) }
								.set { input_kraken2_subsets_tuples }

				input_kraken2_subsets_tuples.view{ "2_subsets_tuples: $it" }

				KRAKEN2( input_kraken2_subsets_tuples )

				KRAKEN2.out.out_list.view { "5_kraken_out: $it" }

				KRAKEN2.out.out_list
					.collectFile(name: "${workflow.runName}.kraken2_out_all.csv", newLine: true)
					.view { "6_after_collectFile: $it" }
					.splitCsv(header: false, sep: ',', strip: true)
					.view { "7_after_splitCsv: $it" }
					.map { row -> tuple(row[0], file(row[1]), file(row[2])) }
					.view { "8_kraken2_out_list_ext: $it" }
					.set { kraken2_out_list_ext }
									//parse the list output by KRAKEN2 to again sampleID and unclas or clas pair
				kraken2_out_list_ext
				    .map { it ->
				        params.USE_KRAKEN_CLASSIFIED == "disable"
				            ? tuple(it[0], it[1])   // unclassified
				            : tuple(it[0], it[2])   // classified
				    }
				    .set { kraken_out1 }
				kraken_out1.view { "9_kraken_out1: $it" }
j
It's likely due to something like your reads coming through the channel in a different order, causing the input to Kraken2 to invalide to cache on resume
You could you probably test this by run once,
then resume
And when it restarts the process, you go into the work directory of both the completed run and the resumed run, and check your input
csv
file
You likely need to collect and sort the lines in the input file
this 2
(note this is a suspicion: I'm not 100 rainbow this is the solution)
n
Also note that even though it is clearly much more appealing than using
=
, usage of
.set{}
is discouraged (https://docs.seqera.io/nextflow/reference/operator#set)
c
Thanks James. That was my suspision too. The input csv files (after subsetted into batches) are identical among runs, so this might not be the reason for it. But i realize that they do have different names, I will set a uniform name and try again I will use = than .set{}. Thanks Niklas
Having the identitcal name solves the resumability issue indeed. 🌟 Thanks for the input ❤️
🎉 1