RADinitio variant in vcf position 0 workaround

[Gwen Ellis]

Hi all,

I am a new RADinitio user trying to use it to complete a power analysis for my anticipated library prep and bioinformatics pipeline based on different demographic scenarios. 

I am using the RADinitio wrapper, but in two simulations it has thrown this error:

[ValueError: A variant position of 0 was found in the VCF output, this is not fully compliant with the VCF spec. If you still wish to write the VCF please use the "allow_position_zero" argument to write_vcf. Alternatively, you can increment all the positions by one using "position_transform = lambda x: 1 + x" or coerce the zero to one with "position_transform = lambda x: np.fmax(1, x)"]

I've looked at the tskit thread for this, but since I'm using the wrapper, I'm unsure how to implement this into my code. Ideally, I would like to either filter out any variants at position 0 (preferred option) or coerce the position from zero to one as each chromosome is simulated. I would be okay with not using the wrapper and using RADinitio with python, but based on the sample code in the radinitio manual, I'm not sure where to implement this since its a tskit error nested in a msprime error nested in a RAdinitio error. 

Any help or clarity on this would be greatly appreciated! Script and output log are attached. I am using python version 3.10, msprime 1.3.4, decoratio 1.0.5, numpy 1.23.5, and scipy 1.10.0

Thank you for the help!

Gwen

[Angel Rivera-Colón]

Greetings Gwen,

Thanks for using RADinitio.

Very interesting. If I am not mistaken, we had a check for this in the RADinitio code, but, as you said, this is coming straight from tskit. I need to dig around a bit more to check where to make the fix. Same with providing a suggestion with the Python interface. At this moment, I just don't know what to change.

In the meantime, since the error happens running the msprime simulation, you could just run RADinitio again using the --make-population flag. This just runs the msprime simulation, meaning you will get new mutations which might not trigger this position 0 error, since the location of the mutations should be effectively random. If you get the same error every time please let me know. Once you finish an msprime simulation successfully, you can continue with RADinitio --make-library-seq, which will simulate the rest of the RAD library. This is a temporary workaround, of course, but it should at least allow you get results while we fix things on our end.

Hope this works in the meantime. Please let me know if you have any questions.

Thanks,

Angel