Large input taking a days to process

[ammaraziz]

Hi again,

I am running MiCall version v7.17.1 through docker, the command I used:

micall folder --project_code HIVB --skip trim.censor -s

The sequencing experiment was hybridisation using the TWIST HIV library run on a P1 NextSeq paired end 150bp. The data generated for most of the samples was heaps (something like 5 million paired end reads on average, unprocessed). These samples with reasonable read numbers finished within hours.

But two samples were nearly triple the number of reads (15 million) are still being processed ~ 4 days later.

The machine I am running MiCall on has plenty of ram (2TB) and cpu (256).

The output (just one sample):

2025-04-24 07:13:42,554[INFO]micall.link_samples(): Pairing files /data/20_XXXXXX_S20_R1_001.fastq.gz and /data/20_XXXXXX_S20_R2_001.fastq.gz.
2025-04-24 07:13:43,460[INFO]micall.drivers.sample.process(): Processing Sample 20_XXXXXX (12 of 24) ('/data/20_XXXXXX_S20_R1_001.fastq.gz').
2025-04-24 08:35:44,843[INFO]micall.drivers.sample.process(): Running fastq_g2p on Sample 20_XXXXXX (12 of 24).
2025-04-24 11:41:52,503[INFO]micall.drivers.sample.run_mapping(): Running prelim_map on Sample 20_XXXXXX (12 of 24).
2025-04-24 19:49:56,932[INFO]micall.drivers.sample.run_mapping(): Running remap on Sample 20_XXXXXX (12 of 24).
2025-04-25 15:48:24,382[INFO]micall.drivers.sample.process_post_assembly(): Running sam2aln on Sample 20_XXXXXX (12 of 24).
2025-04-25 16:47:30,930[INFO]micall.drivers.sample.process_post_assembly(): Running aln2counts on Sample 20_XXXXXX (12 of 24).

To be clear, the pipeline hasn't crashed, it is still processing the data. In the output directory, there is scratch folder containing tmp processing folders. For the above sample number 20 the folder contains 2674 (and still growing) nuc_read_counts*.csv.

The issue is the large number of input reads but what's interesting is that the input is 'only' 3 times larger than files which took ~4 hours to complete.

I could subsample the reads to a reasonable amount but I think the issue still remains that something is triggering MiCall is go super slow. Possibly related to #214 ?

[Donaim]

Hi @ammaraziz

Thank you for reaching out and testing MiCall on your setup.

Your description definitely sounds like something that we have experienced ourselves.
In our case, the bottleneck was IVA (the assembler that we imported from the internet).
When samples sizes are bigger, the slow-down is much bigger, it is non-linear. This part is expected.
But sometimes IVA would take up to 3 weeks for samples that do not differ that much in terms of size.
I think that you are witnessing this second situation.

This hard-to-predict behaviour of IVA distinguishes it from other parts of MiCall, such as the aln2counts.py, which you referenced.
All of the other parts are more-or-less linear in terms of the input size.
We're working on investigating other assemblers to replace IVA, and track that here.

I am attaching a histogram of run-times for a typical run that we do at the BCCFE. It shows that for samples similar in sizes, the run-times are not. (On the plot: x axis is the running time in seconds, y axis is the number of samples falling into a particular run-time bin)

I hope that this answer explains the issue.
As for suggestions, I can only think of trying smaller samples, or waiting for #665 to be resolved.

[ammaraziz]

Thanks for the information.

Is there a way to turn off IVA (denovo) step or is this used for DRM calling?

Also is there a way to predict which samples will cause IVA to take a super long time? Is it related to sample complexity,

[Donaim]

I just realized that IVA is not used in your run. It is only enabled with --denovo flag.
Sorry for my oversight.

If you are not using denovo version, then it is unexpected to see a much longer running time.

We will investigate this.

[ammaraziz]

No problem at all, thank you for helping with this issue. Please let me know what I can do to assist. I'll ask to see if I can securely share data with you for testing purposes.

Does the pipeline output any run metrics I can share with you? It will hopefully now finish after subsetting the super large files.

Purely anecdotal, I think the issue is related to file size, with larger files taking much much longer. The threshold of 'takes a really long time' seems to be around the ~500MB (per read, so ~1G for paired end libs).

[Donaim]

We have multiple debug modes in MiCall, but none of them produce metrics that can be helpful in this case.
I would run micall with a profiler to see which parts of the analysis take the longest.
Maybe if you could do that, then sharing the inputs wouldn't be needed - we would look at the profiler output and see which places need optimizing.
But also the profiler itself slows down the program.

[ammaraziz]

What profile do you recommend?

edit: I think due to the docker setup profiling might be a tad more complicated than I anticipate.

[ammaraziz]

Hi again,

I did some digging into the MiCall code. I think I've identified the where, but not sure on the why.

I mentioned above that thousands of nuc_read_counts*.csv are created and that the pipeline step is aln2counts. Here we can see the code where these files are being generated:

MiCall/micall/core/aln2counts.py

Lines 1720 to 1721 in 09afa22

	file_prefix = os.path.join(os.path.abspath(dirname),
	'nuc_read_counts')

Later on the InsertionWriter Class is called within a for loop:

MiCall/micall/core/aln2counts.py

Line 450 in 09afa22

for i, row in enumerate(aligned_reads):

The nuc_read_counts*.csv files are generated in a reasonable amount of time, ~4 hours between the creation of the first one and the last one. But then nothing happens (no files are created) for hours (currently ~10 hours). Something downstream from here is taking a very long time.

Also this files are empty:

.

I should probably run a profiler haha.

[ammaraziz]

So I found py-spy can inspecting running python apps given a PID. It looks like big_counter.py is the culprit.

I think (but not sure) those are line numbers at the end of each flamegraph bar.

profile.txt - download, rename to svg then open in your browser, github borks svg files for some odd reason.

Reading the flamegraph, it seems to flicker between these lines:

MiCall/micall/utils/big_counter.py

Lines 16 to 22 in 09afa22

	with open(cache_file_name) as f:
	f.seek(pos)
	line = f.readline()
	if not line:
	break
	pos = f.tell()
	yield line

I also ran strace -p $PID, the system calls to open nuc_read_counts*.csv files, over and over again.

If I understand the code, it is opening a nuc_read_counts*.csv file, reading a specific line, then closing and moving to the next line.

The above points to IO bound operation, especially when the number of files is ~1400. What would the effect be of increasing max_size=5000 to 20k or bigger?

[Donaim]

Thank you for the investigation!
It seems like py-spy is the perfect profiler for this task.
And the BigCounter class does look like a place that can be causing the slow down.
Looking at the code, my hypothesis is that it could be hitting the in-memory cache limit (more than 5000 entries) and then starting to use the disk mostly, when we get to the number of reads that are in the 15 million range.
I will see how the BigCounter can be further optimized, but increasing the max_size to 20k should help. Given that you have 256 GB of RAM, I would put it even higher, maybe 100k?

[ammaraziz]

Given that you have 256 GB of RAM, I would put it even higher, maybe 100k?

The server I am running things on has a ton of ram (2TB), so I don't think it's hitting any limit but instead the hard cutoff is to restrictive.

How can I run MiCall with modified settings to affect the max_size, eg does MiCall take any config file? I'm trying to avoid forking MiCall and having to create a new docker (done this before, makes things difficult down the track when updates get pushed upstream).

[Donaim]

The max_size cannot be changed from a config, it is hard coded in the source code.
So I can only see one way to adjust it - to rebuild MiCall docker image.
In the meantime I will make sure to investigate this issue and optimize the BigCounter

If you wish to rebuild the image yourself, this should be easy:

git clone https://github.com/cfe-lab/MiCall
cd MiCall
git checkout v7.17.1
docker build --tag micall:adjust_max_size .

then use the micall:adjust_max_size name when running it.

[Donaim]

Hi @ammaraziz

We have run some MiCall experiments with various input sizes, but weren't able to reproduce the exact behaviour that you reported.
For us, the run time scaled proportionally to the size of the input.
However, the error rate of the input reads matters a lot - higher error rate corresponds to much longer processing times.

Have you been able to deal with the problem on your end?

[ammaraziz]

Hi @Donaim,

Apologies for the slow response.

The solution I came up with was to randomly subsample to a reasonable read count per sample. This isn't ideal, because I am effectively throwing out data, but it's the solution that worked.

I think exposing max_size value to the cli would help to reduce time, the default is very strict and could be expanded.