gaynorr / AlphaSimR

R package for breeding program simulations
https://gaynorr.github.io/AlphaSimR/
Other
42 stars 18 forks source link

Issues with IBD tracking #114

Closed janaobsteter closed 1 year ago

janaobsteter commented 1 year ago

Describe the bug We've noticed something weird going on with IBD tracking when working with drones (double haploids!) in SIMplyBee. Since they are double haploids, all of their offspring should receive the exact same genetic material. This is true when looking at the IBS, but not when looking at the IBD! The IBD relationship with the father is hence everything between 0 and 0.5, where in fact it should be 0.5 for everyone. I've managed to reproduce the bug in pure AlphaSimR code, although it behaves a bit differently. I've created a single mother and a double haploid father and then created offspring with makeCross2. Within a single cross, all of the offspring inherit the same IBD chromosome from the father, but it's not always the same! I am hence running 100x crosses to show the issue.

Steps To Reproduce

# AlphaSimR IBD tracking problem
founderPop = quickHaplo(nInd=10, nChr=1, segSites=1000)
#Set simulation parameters
SP = SimParam$new(founderPop)
# Track the pedigree
SP$setTrackPed(TRUE)
# Track the recombination
SP$setTrackRec(TRUE)

#Create population
basePop = newPop(founderPop, simParam=SP)

#The first basePop individual is the mother
mother <- basePop[1]
# The second basePop individual will serve as a "mother" to create DH fathers
father <- makeDH(pop = basePop[2], nDH = 1, keepParents = FALSE, simParam = SP)
# Check whether it's truly a DH
fatherIBDHaplo <- pullIbdHaplo(father)
table(fatherIBDHaplo[1,] == fatherIBDHaplo[2,])
# Is it a true double haploid!
table(fatherIBDHaplo[1,]) # This is father's IBD haplotype - should be passed "unrecombined" to the offsprign

# Make some offspring with randCross2 (repeat it 100X times to catch the error)
i = 1
haplo <- vector(mode = "list", length = 100)
while (i < 100) {
  offspring <- randCross2(females = mother, males = father,
                          nCrosses = 10, nProgeny = 1, balance = FALSE)

  # Check how the offspring matches the DH father
  for (id in offspring@id) {
    haplo[[i]] <- list(table(pullIbdHaplo(offspring)[paste0(offspring@id[1], "_2"),]))
  }
  i = i + 1
}
# 3 and 4 are IBD haplotypes inherited from the father - should be all the same!
table(unlist(sapply(haplo, FUN = function(x) x[[1]]["3"])))
table(unlist(sapply(haplo, FUN = function(x) x[[1]]["4"])))

Expected behavior All of the offspring should receive the same haplotype from the DH father, but they don't. I do not know whether the issue is in the DH or in tracking in general, since I can't check whether the tracking from the mother to offspring is correct, since it recombines and it's hard to catch. However, mother-offspring IBD relationships in SIMplyBee seem to be fine (0.5).

gaynorr commented 1 year ago

Thanks for sharing this. I'll take a look at it right away.

gaynorr commented 1 year ago

This issue should now be fixed in the devel branch.