Relatedness estimates - packages suggestions & "related" format obtained through gl2relate

[Gabriella Scatà]

Hi everyone,

this is not a question about dartR, apart from the last bit on the related format obtained through the gl2relate function.

However, it might be related as there are several functions in dartR allowing conversions to packages computing relatedness (i.e. gl2relate, gl2demerelate, gl2gds), but I am finding that most of these packages have functions that don't work anymore or at least not with large SNPs datasets.

Which R package do you use to estimate relatedness, especially if you have small population sample sizes?

So far I tried:

- related --> RStudio crashes every time I try the coancestry() function - regardless of the size of the SNPs dataset (even if i use a 50SNPs subset)

- Demerelate --> works only with 100 SNPs  - any dataset with more loci it gives an error, and the function loci.test() in Demerelate does not work either with any dataset size

- SNPrelate is the only one which seems to work so far, but I still have to read the relevant literature as the method it relies on might not be accurate enough for small sample sizes (population n~ or < 20 individuals)

I would like to try the COANCESTRY software which has a bias correction for small populations, but the "related" file format i obtain from the gl2relate function, which is seen as a data frame by R, does not seem to be written into a normal .csv file when I use write.csv / write_csv / write.table...is there a reason for it?

The software COANCESTRY needs a pure text file format, using either comma, tab, or white space as delimitators, so I am assuming a normal .csv file or table file would be ok. But, the software does not seem to accept the related format file written into a .csv file.

Any suggestion or info any of you has on this would be highly appreciated.

Thanks everyone for your input and help as usual!

Cheers,

Gabriella

[Jose Luis Mijangos]

Hi Gabriella,

Please find attached a document with some code to run four kinship approaches:

1. dartR function gl.grm.network (read function documentation).

2. Function to convert genlight to colony format written by Sunnuck’s group (see attached article).

3. Code to convert genlight to sequoia R package format (see Huisman, Jisca. "Pedigree reconstruction from SNP data: parentage assignment, sibship clustering and beyond." Molecular ecology resources 17.5 (2017): 1009-1024.)

4. Code to convert genlight to related R package format which runs the approach of the coancestry program (see Wang, Jinliang. "COANCESTRY: a program for simulating, estimating and analysing relatedness and inbreeding coefficients." Molecular ecology resources 11.1 (2011): 141-145.)

Cheers,

Luis 

[Gabriella Scatà]

Hi Luis,

thank you so much for taking the time to outline a few approaches to estimate relatedness. This is so helpful!

I will go through them slowly and let you know if I have questions/issues.

The "related" package is the only one from the ones you listed that I tried, and unfortunately R crashes whenever I run the "coancestry()" function of the "related" package, no matter the size of the dataset (25/50/more SNPs).

I will test the other approaches!

Thanks a lot,

Best,

Gabriella

--
You received this message because you are subscribed to a topic in the Google Groups "dartR" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dartr/YBgRrEqMzgM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dartr+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dartr/f8df2596-bb77-4c4a-b01f-50db69eca873n%40googlegroups.com.

[Jose Luis Mijangos]

Hi Gabriella,

You can try following the instructions in the link below to install the related package.

https://github.com/timothyfrasier/related 

Cheers,

Luis 

[Gabriella Scatà]

Hi everyone,

I was just wondering whether there is any function in dartR to convert to .bam files?

I would like to use the IBSRelate method as I have population structure in my dataset (in the form of hybridising species) and I would like to use a relatedness approach that can be robust to population stratification with geographic population sample sizes around 20 individuals each.

I though the method developed by Waples (2018) --> https://onlinelibrary.wiley.com/doi/10.1111/mec.14954 would be ideal for my case.

Do you have any suggestion on this matter in case?

Thank you!

Best,

Gabriella

[Jose Luis Mijangos]

Hi Gabriella,

We don't have a function to convert to bam files. 

You could try a new method developed by Wang that is robust to population structure and population size:

Wang, J. (2022). A joint likelihood estimator of relatedness and allele frequencies from a small sample of individuals. Methods in Ecology and Evolution, 13(11), 2443-2462.

library(devtools)
install_github("green-striped-gecko/dartR@dev")
library(dartR)
# Download the program from here:
# https://www.zsl.org/about-zsl/resources/software/emibd9
# For Windows, install the program and then move the following files to your
# working directory: "EM_IBD_P.exe", "impi.dll" and "libiomp5md.dll".
# For Mac move the file "EM_IBD_P" to the working directory.

res <- gl.run.EMIBD9(platypus.gl)

Cheers,

Luis 

[Gabriella Scatà]

Hi Luis,

and thank you so much for this.

I had actually seen one of your previous responses to a similar post, and had already tried to download and install the EMIBD9 software.

To obtain a Genotype data, I had followed this steps:

<

m = as.matrix(my.genlight)   # obtain the matrix with genotype data coded as 0/1/2 with nInd x NLoc (rows x col) = DataFormat2 of EMIBD9 software

m[is.na(m)] = 3   # convert all NA (missing values) to "3" (as required from the data format of EMIBD9)

m= as.data.frame(m) # convert to dataframe in order to export the file to .csv to load it in EMIBD9

>

However, when I load the file in EMIBD9 and set up all the parameters (es. n individuals, Inbreeding=1, etc...), and I click "save", the program tells me "Error: expected 425 individuals, but only 1 detected".

I tried to remove the colnames of the dataframe, since for the similar software COANCESTRY the colnames with all the names of the loci had to be removed from the input file.

But that did not make any difference for the EMIBD9 software.

I then tried running the EMIBD9 software from R with the approach you suggested (I moved EMIBD9 files to the working directory):

- install dev@dartR version (current as of 15 nov 2023)

- run res <- gl.run.EMIBD9(my.genlight)

But unfortunately, the institution I work for has the Firewall/IT Security system blocking any external program run through R...I do have administrators rights so I can run EMIBD9 independently outside of R, but I am not sure how to overcome this issue from within R.

If you have any input, I would really appreciate it.

Otherwise, if you have any suggestions on why when I run the EMIBD9 software externally from R, I get the error message reported above, I would also appreciate any feedback.

Thanks again for all your help,

Best,

Gabriella

[Jose Luis Mijangos]

Hi Gabriella,

Please find attached a function to convert a genlight object to EMIBD format. Please download it, unzip it, and move it to your working directory. Then try the code below. 

library(dartR)
source("gl2EMIBD9.R")
gl2EMIBD9(platypus.gl)

Cheers,

Luis 

[Gabriella Scatà]

Hi Luis,

thank you for the alternative solution.

I tried and followed your instructions:

- downloaded the function you sent, extracted and moved it to my wd

- run the code source("gl2EMIBD9.R") (I can see it now in my list of functions)

- loaded dartr (I tried with both dartR versions: loading either the current dartR version or the dev version 2.9.7)

- run gl2EMIBD9(my.genlight) either like this or saving the results into an object (res = gl2EMIBD9(my.genlight)

But nothing seems to happen.

When I save the res into an object and I try to view it, it turns out as "type" = NULL, Value = pairlist of length 0.

Any idea?

Today I will try to reach out to our IT for alternative solutions.

Thanks again!

Best,

Gabriella

[Jose Luis Mijangos]

Hi Gabriella,

The two input files should have been saved in your working directory.

Cheers,

Luis 

[Gabriella Scatà]

Ooh, yes I see them now!

Sorry, I hadn't even thought of checking there!

Thanks a lot!

I'll test it in the EMIBD9 software now!

Best,

Gabriella

[Gabriella Scatà]

Hi Luis,

I managed to run the analysis in EMIBD9.

I just got the results, however, they are illegible.

I know I can open them in EMIBD9 software to visualize them, but I want to run a mantle test with these and other kinship estimates in R and the results are not even in a table format...the kinship values are just listed one after the other in a row (similar to the input genotype file).

Do you have by chance also a function that converts/extracts the results from the EMIBD9 output file?

Thanks again!

Best,

Gabriella

[Jose Luis Mijangos]

Hi Gabriella,

Download the folder, unzip it, move files to the working directory, and try:

library(dartR)
source("gl2EMIBD9.R")
source("gl.read.EMIBD9.R")
gl2EMIBD9(platypus.gl)
res <- gl.read.EMIBD9(file.name="EMIBD9_Res.ibd9")

Cheers,

Luis 

[Gabriella Scatà]

Hi Luis,

thank you so much! You guys are my saviours!

I managed to open the .ibd9 results file from EMIBD9 software through the function you sent (gl.read.EMIBD9).

However, the sample IDs names are not there anymore (I think they were not even supposed to be there in the input file for the EMIBD9 software), and also the order of the samples in the matrix obtained through gl.read.EMIBD9 is all over the place (i.e. 1, 10, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 11, 110, 111,...etc...)...

So, it's hard to know i) if the numbers 1, 2, 3 etc...correspond to the original order of samples in my genlight object es.: 18shk0006, 18shk0007, etc...) and ii) as the order of samples is not as the original one, it will be difficult to compare these results with those obtained for example from KING robust with a mantel test as the order of the samples should be the same in the matrices being compared.

So, my questions are:

1. does the numbering 1,2, 3, etc...to 425 correspond to the order of my original samples in the genlight object?

2. If so, then I can find a way to re-order the samples in the matrix obtained with gl.read.EMIBD9 so that i can obtain the original samples order in my genlight object?

Thank you again for your constant assistance and feedback!

Best,

Gabriella

[Bernd.Gruber]

Hi Gabriela,

 

We are about to have a new function implemented in the new package dartR.captive. It is not 100% ready yet and yes the current function renames the names in the emibd9 file to avoid problems with spaces in indnames.

 

The new function will return a matrix with the same indNames as the genlight object (and in the same order as provided), but as mentioned will take another week or so until tested.

 

I will let you know when finished.

 

Cheers, Bernd

--
You received this message because you are subscribed to the Google Groups "dartR" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dartr+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dartr/7d6c8b64-b8cd-4d97-b01d-61ffd4d6c61dn%40googlegroups.com.

[Jose Luis Mijangos]

Hi Gabriella,

 

I suspect that you used the first function I sent (gl2EMIBD9). I fixed that names issue in the second function I sent you.  

Cheers,

Luis 

[Gabriella Scatà]

Ahhh ... ok that's why you sent the gl2EMIBD9 function again...so it means I have to re-run the input file obtained with the new function you sent me (gl2EMIBD9.R) in the EMIBD9 software to obtain a new output file with the correct order of samples?

And, just to be clear, the order of samples - that is the numbering 1-425 -  in the .ibd9 output file, will then correspond to my original samples order in the genlight object?

Thanks again so much to be so on top of things!!

Best,

Gabriella

[Gabriella Scatà]

Hi Bernd,

sounds great :)

Looking forward to test it as well.

Thank you to all of you for your constant assistance and great work!

Best,

Gabriella

[Gabriella Scatà]

Hi Luis, 

just to let you know everything worked great! 

Thank you so much for sending those functions over!!

Best,

Gabriella

[Manuela Cascini]

Hi Luis,

I tested EMIBD9 on one of my dataset using the function  gl.run.EMIBD9 and had no issues. 

However, I am confused by the results. Specifically, I anticipated that self-comparisons would yield consistent values across the dataset. Instead, the values range from 0.59 to 0.82.

Have you came across similar results? Do you know the reason for this?

Best,

Manuela

[Jose Luis Mijangos]

Hi Manuela,

The diagonal of the relatedness matrix is the self coancestry. You can find more information in the first paragraph of page 4 (highlighted) in the attached article. 

Cheers,

Luis

[Elie McCarthy]

Hi Luis,

Is there any way to generate confidence intervals for the kinship values using EMIBD9 like there was for coancestry?

Many thanks for your help,

Elie

[Jose Luis Mijangos]

Hi Elle,

Unfortunately, there isn’t a built-in option to calculate confidence intervals for EMIBD9. However, both Milligan’s method (dyadml) and Wang’s triadic method (trioml) in the Coancestry package are highly accurate and do allow for confidence interval estimation.

Cheers,

Luis