RELERNN SIMULATE issue with vcf?

[dylandebaun]

RELERNN SIMULATE issue with vcf?

[jradrion]

@andrewkern I should have some time to take a look at this over the weekend. Hopefully there will be an easy fix.

[andrewkern]

hi there-- can you share the command line you used to run the program?

[dylandebaun]

Oh yes, here it is. I've also run it with --forceDiploid I believe and gotten the same error.

DIR="$PWD"
MU="2.4e-9"
ReLERNN_SIMULATE --vcf Leioheterodon_madagascarensis_B_biallelic.vcf --genome Leioheterodon_madagascarensis_B_biallelic.vcf.bed --phased --projectDir ${DIR} --assumedMu ${MU} --upperRhoThetaRatio 5 --assumedGenTime 2 --nTrain 12800 --nVali 2000 --nTest 100 --nCPU 40

[andrewkern]

A quick guess @jradrion, given that this is from simulate. The msprime API moved such that by default sample number is diploid individuals and not haploid chrom number-- I wonder if that's messing stuff up?

[jradrion]

@andrewkern This isn't an issue with msprime output as the error occurs prior to simulation. Looks like there may be an issue with differing ploidy levels across the chromosomes (scaffolds/contigs/etc). The code creates a mask from missing genotypes for each chromosome and then concatenates the masks using np.concatenate(). The error is being cause because the shapes on the arrays being concatenated differ (some arrays are created using 2n=6 but at least one chromosome appears to be being interpreted as 1n=3). The error message is calling out the array at index 6, which from the stdout I think will correspond to the following chromosome, 8739_RagTag. The genotypes for that chromosome aren't in the first 50k lines of the attached .vcf file. @dylandebaun can you provide the .vcf data for it? I'm traveling the next three days and will not have a chance to look at this, but I'll try to get back to you when I have the chance. You could also try excluding this chromosome and see if you get the same error. Right now I'm not sure where that differing number of haplotypes is being introduced.

[dylandebaun]

Ah that's what the array was referring to! Yes, that was the issue, one of my vcf filtering steps removed that chromosome from the vcf file, which might explain why it wasn't caught in the original ReLERNN error that told telling me which of my chromosomes were hemizygous. I reran it without that chromosome in the bed file and again with the correct vcf and both ran to completion.
Thank you so much for the help!

[andrewkern]

good detective work @jradrion. I'm going to close this issue @dylandebaun.