Assistance with filtering and relatatedness analysis

[Molly Barlow]

Hi all, I'm attempting some simple relatedness and inbreeding analysis of a fenced population of marsupial, but my relatedness networks are indicating a much higher degree of relatedness than they should be (i.e. suggesting they're all identical). I was also hoping for some advice on my filtering procedure. 

I've sent my code directly to Luis. 

Cheers,

Molly 

[Jose Luis Mijangos]

Hi Molly,

Thanks for sending your code, your filtering already looks great. Just a few comments.

- Remove individuals with low call rate first, with gl.filter.callrate(x, method = "ind"). Doing this first keeps SNPs you would otherwise lose.
- Use gl.filter.rdepth(x, lower = 5, upper = 1000). Loci with low depth are often miscalled as homozygotes, so this gives you more accurate heterozygote calls.
- Use gl.filter.maf(x, threshold = 0.05). Rare alleles add noise to relatedness estimates and are frequently sequencing errors.

On the relatedness estimates themselves, the thing to keep in mind is that relatedness is relative, not absolute. It is estimated against a reference population, which in practice is your own dataset, and estimators assume that reference population is outbred and randomly mating. If several samples in the dataset are related or inbred, estimates are biased upwards. If the dataset has samples from two differentiated populations, they are biased downwards.

Methods handle this differently. gl.grm builds a genomic relationship matrix (GRM) centred on your sample, so 0 is not an unrelated pair but a pair of average relatedness for your dataset. In a dataset like yours, descended from a few founders, the average pair already shares a lot of ancestry, so a GRM value of 0 already sits above true unrelatedness.

I would use EMIBD9 instead, via gl.run.EMIBD9 with Inbreed = TRUE. It estimates allele frequencies and relatedness coefficients jointly by maximum likelihood, which partly corrects for the problem above.

Also note that relatedness and kinship have different definitions. Relatedness is the expected proportion of alleles that are identical by descent between two individuals, so full siblings have a theoretical relatedness value of 0.5. Kinship is the probability that two alleles, one drawn at random from each individual, are identical by descent, so full siblings have a theoretical kinship value of 0.25 (i.e. kinship is half of relatedness).

 EMIBD9 reports kinship, and you can feed it straight into gl.grm.network, as shown below.

library(dartRverse)
library(viridis)
programs_path <- "path_to_folder_containing_EMIBD9"

t1 <- platypus.gl
EMIBD9 <- gl.run.EMIBD9(t1,
                        emibd9.path = programs_path,
                        Inbreed = TRUE)
grm_sib <- gl.grm.network(G = EMIBD9$rel,
                          x = t1,
                          categorise = TRUE,
                          node.size = 8,
                          palette_discrete = viridis::viridis_pal(),
                          kinship.threshold = 0.125,
                          legend.title = "Sampling locations",
                          title = "")

You can download the EMIBD9 binary from here: https://github.com/green-striped-gecko/dartRverse/tree/main/binaries

Happy to help further if anything is unclear.

Cheers,
Luis