I have two questions about the PopGenome package for R.
- When using the
readVCFfunction, it is possible to read several chromosomes at the same time? If not, it is possible to merge several
GENOME.classobjects? (my VCF files are "small", I don't have memory limitations).
- Coalescent simulations are computed correctly well when applied to a whole chromosome but fail when applied to a subset of sites: i.e. exons, introns, sliding window analysis... I copy the code for clarity:
For a certain chromosome (works well!):
GENOME.classChr1 <- neutrality.stats(GENOME.classChr1) MS.classChr1 <- MS(GENOME.classChr1,thetaID="Tajima",neutrality=T)
For exons (fails!):
GENOME.classChr1.splitEx <- neutrality.stats(GENOME.classChr1.splitEx) #this command works well MS.classChr1.splitEx <- MS(GENOME.classChr1.splitEx,thetaID="Tajima",neutrality=T) #ERROR: Error in 0:nsamtot : argument of length 0
How can I solve these issues?
Thanks for your suggestion about reading different chromosomes at the same time. In fact is a good idea and should work, the problem is that my VCF files are huge and R session is aborted when reading data in this way...
I've solved the issue of the coalescent simulations by running first the ms program to get the simulations (which allows you to modify default parameters) and read the output file in popgenome using the
readMS()function. Once you have read the simulations file you can compute different statistics, for instance using
Hope it helps!