Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Kindepth for double marriage #33

Open
LouisLeNezet opened this issue Nov 14, 2023 · 3 comments
Open

Kindepth for double marriage #33

LouisLeNezet opened this issue Nov 14, 2023 · 3 comments

Comments

@LouisLeNezet
Copy link
Contributor

Hi @sinnweja,

I'm nearly done for the new package.
I've recently discovered a subtle use case where the kindepth is not nicely computed:

Here is the code:

df <- data.frame(
    id = 1:12,
    dadid = c(0, 0, 1, 0, 0, 0, 3, 3, 5, 5, 7, 10),
    momid = c(0, 0, 2, 0, 0, 0, 4, 4, 6, 6, 9, 8),
    sex = c(1, 2, 1, 2, 1, 2, 1, 2, 2, 1, 1, 2)
)
pedi <- with(df, pedigree(id, dadid, momid, sex))
kindepth(pedi, align = TRUE)
# 0 0 1 1 0 0 2 2 1 1 3 3

We get
image

What would be better is to to have a kindepth of

kindepth(pedi, align = TRUE)
# 0 0 1 1 1 1 2 2 2 2 3 3

plot

The problem seems to come in the loop that update the agood ancestors.
In the case of a double marriage like in the previous pedigree, the loop will iterate and will add the individuals as follow:

# Initial depth
0 0 1 0 0 0 2 2 1 1 3 3

# id 4 single marryin corrected
0 0 1 1 0 0 2 2 1 1 3 3

# Selected pairs is 9, 7
# abad is
9 5 6
# agood is
7 1 2 3 4

# Then by iteration agood become
# + spouse
7 1 2 3 4
# + ancestors
7 1 2 3 4
# + kids
7 1 2 3 4 8
# + spouse
7  1  2  3  4  8 10
# + ancestors
7  1  2  3  4  8 10  5  6
# + kids
7  1  2  3  4  8 10  5  6  9

The individuals 10, 5, 6 and 9 should not be added to be able to modify them.
What would be the best way to control for this case ?

Best,
Louis

@LouisLeNezet
Copy link
Contributor Author

What I've found to control that is to do the following:

if (length(abad) == 1 && sum(c(dads, moms) == bad) == 1) {
...
else {
    ## Ancestors of the 'good' side
    agood <- c(good, ancestors(good, midx, didx))
    ## If individual already in agood they aren't really bad
    abad1 <- abad[!abad %in% agood]
    ## For spouse chasing, I need to exclude the given pair
    tdad <- dads[-who]
    tmom <- moms[-who]
    ## Get all individuals affiliated to agood
    while (TRUE) { ... }
    ## Update agood but only if not in abad1
    agood <- agood[!agood %in% abad1]
    ...
}

By doing so we keep a good kindepth measure in all case I have for now.
What do you think ?

@sinnweja
Copy link
Member

sinnweja commented Nov 14, 2023 via email

@LouisLeNezet
Copy link
Contributor Author

Hi Jason,

I've added this change to the current "Pedixplorer" package for Bioconductor.
But I can do a PR to kinship2 with this change, if you want.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants