bcftools mpileup doesn't open the BAM index file when used with -b and large number of samples

[arturotorreso]

Hi,

I am trying to call variants in a list of BAM files that contains 2600 bam files from bacteria.

The command is as follows: bcftools mpileup -r $region -a AD,INFO/AD,ADF,INFO/ADF,ADR,INFO/ADR,DP,SP -pm3 -b $bamList -E -f $REF -Q 0 -Ou | bcftools call -mv -f GQ,GP -p 0.99 -Oz -o $outputVCF

Bcftools mpileup gives me the following error:

[mpileup] fail to load index

The index is indeed where it should be, and it's not corrupt. The file that failed is in the line 2044 in the list of bam files. I shortened the list of bam files to 500, containing the file that previously failed, and then it ran perfectly. I shuffle the files in the long list, and bcftools mpileup failed in a different file, for the same reason, also in the line 2044 of the bam list. This indicates me that the files themselves are not the problem.

Is there a limit on how many samples should be run simultaneously? Is it given the wrong error for another problem?

Thank you for your help, Arturo

[pd3]

Try to run ulimit -n to see if there is a system limit on the number of simultaneously open files. The limit can be changed by running eg ulimit -n 9999.

If that's not possible because of permissions (bash: ulimit: open files: cannot modify limit: Operation not permitted), edit the following files and add these lines:

sudo gedit /etc/pam.d/common-session
        session required        pam_limits.so

sudo gedit /etc/security/limits.conf
        *      -       nofile          9999

Then log out and then log in back. The command ulimit -n should then print the value 9999.

This is distribution-specific, but should give you enough hints to solve on your own.

[arturotorreso]

Hi Petr,

that was indeed the problem! I didn't think of it.

I changed the ulimit -n and now it runs perfectly.

[Elakhdar]

Hi Petr, I still have a similar issue however the number of simultaneously open files was changed to 9999

!!ERROR!! bcftools mpileup -a AD,ADF,ADR -B -q 40 -Q 18 -C 50 -O u -r chr03 -f result/10_ref/IRGSP-1.0_genome.fasta result/20_bam/cultivar.filt.bam result/20_bam/bulk.filt.bam | bcftools call -vm -f GQ,GP -O u | bcftools filter -i "INFO/MQ>=40" -O z -o result/30_vcf/mutmap.chr03.vcf.gz >> result/log/bcftools.chr03.log 2>&1

[pd3]

@Elakhdar please open a new issue, doubt this is the same problem since you work with one input file, not thousands of files. You don't show the exact error message so can't be of much help here. Make sure you try with the latest version of bcftools+htslib, oftentimes it gives more informative error messages than older versions.