suggestion: Consider re-organizing this block to something like:

if (params.stat_mode == 'library') {
    filename_id = lib_id
    outdir = ...
    log_suffix = ...
} else if (params.stat_mode == 'readgroup') {
    filename_id = rg_id
    outdir = ...
    log_suffix = ...
} else {
    filename_id = sm_id
    outdir = ...
    log_suffix = ...
}

output_filename = generate_standard_filename("SAMtools-${params.samtools_version}",
            params.dataset_id,
            filename_id,
            [:])

and then using addParams when importing the process to set stat_mode.

This reduces the amount of duplicated code (ex. the standard filename generation function is only coded once) and also removes the need for the process itself to have to be imported with a specific name and instead be controlled by a parameter set during import.

I forgot to change this. In a separate PR I will make this a list of readgroups and add the -r within the processes

Some of the GIAB samples for example have > 1000 readgroups. I think I should impose a maximum number of readgroups, and I suppose libraries, that will be processed separately by SAMtools stats. But for FastQC they will all be processed separately as this is equivalent to running FastQC on the FASTQ files.

I added parameters for maximum number of readgroups and libraries to process separately with defaults = 20.

I also parameterized a choice of levels to run FastQC.

@sorelfitzgibbon
question: I think I understand what this command tries to archive but can you elaborate on this a bit (can this -F 0x900 be optional, etc)? Also, have you tried to incorporate threading? This would be probably really slow for large samples. Is it possible to use fastqc -format bam with threading for the sample level run?

Thanks @tyamaguchi-ucla

• -F 0x900 removes secondary and supplementary alignments. It is actually the default but I thought I'd leave it here so there is no question that its being done.
• Surprisingly fastqc -format bam does not allow per readgroup analysis.
• I have not tried to incorporate threading. All three commands allow threading. I will try it and test speedup.

UPDATE 2: After discussion with @yashpatel6:

• Since, within metapipeline, this pipeline will eventually be used to gate the running of further pipelines, realtime will be more important than cpu hours, thus the default will be 2 cpus. More than 2 cpus provides very little time advantage for the large (more typical size) readgroup.
• cpus will be adjusted down to 1 for input with large numbers of readgroups.

UDPATE: I now think it would be better to avoid threading entirely, at least when running in metapipeline.

• Testing on a sample with 112 readgroups took twice as long with 2 threads vs 1 (9.5 hours vs 4.7 hours).
  /hot/software/pipeline/pipeline-generate-SQC-BAM/Nextflow/development/unreleased/sfitz-add-mosdepth/slurm-184021.out
/hot/software/pipeline/pipeline-generate-SQC-BAM/Nextflow/development/unreleased/sfitz-add-mosdepth/slurm-184029.out
• The per cpu benefit of 2 threads on samples with fewer readgroups is minimal, if even real.

Below are the results for running with various threads on SGSALASC000015-T011-O04-F.bam (185 GB).

• This sample has five readgroups with one being much larger than the others.
• All output here: /hot/software/pipeline/pipeline-generate-SQC-BAM/Nextflow/development/unreleased/sfitz-by-readgroup
• It looks like the sweet spot is 2 threads

Speedup (single thread realtime/realtime)

1 thread; full run CPU HOURS = 9.3

2 threads; full run cpu hours = 9.2

3 threads; full run cpu hours = 12.2

4 threads; full run cpu hours = 12.4

6 threads; full run CPU HOURS = 16.7

12 threads; full run CPU HOURS = 37.0

-F 0x900 removes secondary and supplementary alignments. It is actually the default but I thought I'd leave it here so there is no question that its being done.

Thanks for the benchmarking results. I wasn't aware of this default. Do we want to include them though (i.e. -f?) to get a complete result for the RG level FASTQs, etc?

No, I don't think we do want to include them. They are duplicates of reads that are already represented. Removing them creates the output that directly matches the original FASTQs.

@sorelfitzgibbon (suggestion/question: non-blocking) You've generated lots of benchmarking results and we should post them on GitHub discussions as well. Also, do you typically adjust for caching performance when conducting benchmarking?

do you typically adjust for caching performance when conducting benchmarking

Is there a recommended way of doing that? I ran almost all of it on nodes that were already running and had space available. They were run at more or less the same time. I could run them all concurrently on a single F72, although I assume even that background will change as jobs start to finish.

Caching adjustments are hard to control in general; a way around it could to be make the tests totally independent of the cache by copying the input files to /scratch and then run all benchmarking runs using files under scratch.

I will do some cache independent benchmarking.

Caching adjustments are hard to control in general; a way around it could to be make the tests totally independent of the cache by copying the input files to /scratch and then run all benchmarking runs using files under scratch.

Yeah, I was going suggest this. Maybe worth creating a confluence page for general benchmark recommendations? We've discussed pre-caching, etc with OHIA and MSFT a while ago but using /scratch would be the easiest solution.

I'm ok with maintaining the three sections and using run_statsReadgroups_SAMtools

suggestion: Do we also want to add the max libs parameter here?

@yashpatel6 I left it off in the interest of streamlining the template since library numbers that are large enough to be an issue are currently unlikely. It is in the README. Do you want me to add it?

For consistency, I would suggest adding it