GRAPE: genomic relatedness detection pipeline

Classifying the degree of relatedness between pairs of individuals has both scientific and commercial applications. As an example, genome-wide association studies (GWAS) may suffer from high rates of false positive results due to unrecognized population structure. This problem becomes especially relevant with recent increases in large-cohort studies. Accurate relationship classification is also required for genetic linkage analysis to identify disease-associated loci. Additionally, DNA relatives matching service is one of the leading drivers for the direct-to-consumer genetic testing market. Despite the availability of scientific and research information on the methods for determining kinship and the accessibility of relevant tools, the assembly of the pipeline, which stably operates on a real-world genotypic data, requires significant research and development resources. Currently, there is no open source end-to-end solution for relatedness detection in genomic data, that is fast, reliable and accurate for both close and distant degrees of kinship, combines all the necessary processing steps to work on a real data, and is ready for production integration. To address this, we developed GRAPE: Genomic RelAtedness detection PipelinE. It combines data preprocessing, identity-by-descent (IBD) segments detection, and accurate relationship estimation. The project uses software development best practices, as well as Global Alliance for Genomics and Health (GA4GH) standards and tools. Pipeline efficiency is demonstrated on both simulated and real-world datasets. GRAPE is available from: https://github.com/genxnetwork/grape.


Introduction
Distant relationship estimation has both scientific and commercial applications. Scientific applications may include the identification of monogenic (single gene) Mendelian diseases. 1,2 Also, relatedness detection can be used during data quality control for genome-wide association studies (GWAS), since close relatives should be excluded to ensure that no pair of individuals are more closely related than second-degree relatives. 3,4 Otherwise, GWAS may suffer from high rates of false positive results. The potential commercial application of relationship estimation is utilized by direct-to-consumer genetic testing companies to find possible distant relatives for their customers. 5 The applicability of existing tools and pipelines to both applications is limited. Currently, there is no open-source end-to-end solution for relatedness detection in genomic data, that: (a) is reliable and accurate for both close and distant degrees; (b) includes all necessary processing steps to work with real data; (c) is flexible enough to be adapted for various applications; (d) is user-friendly and ready for production integration. Specifically, for commercial usage, the pipeline should deal with data heterogeneity and efficiently process newly added data samples. Samples can be genotyped with different chips and may have been passed through different quality controls. Typical databases of genetic testing companies contain a lot of data from previously used chips, which need to be combined during the processing. The database may also contain heterogeneous data. For example, UK Biobank dataset 6 contains samples which were genotyped using two different chips.
Driven by the idea of open and user-led innovation, we have developed GRAPE (Genomic RelAtedness detection PipelinE), the first open-source end-to-end solution for relatedness detection that can successfully address the abovementioned difficulties. As a preliminary step, we comprehensively studied available approaches and software instruments for relatedness estimation from genotype data, such as identity-by-descent (IBD) segments detection tools: GERMLINE, 7 IBIS, 8 RaPID, 9 PhasedIBD 10 ; relationship inference tools: DRUID, 11 ERSA, 12,13 KING 14 ; other tools, which may be required during data preprocessing steps, like algorithms of phasing (Eagle 15 ) and genotype imputation (Minimac 16 ). Then, we selected several perspective tools and joined them together into the user-friendly GRAPE pipeline.
GRAPE adapts the best practices for software development, including the Snakemake 17 workflow management system, Conda 18 virtual environments, Docker 19 containerisation, and Funnel task execution service, which implements GA4GH Task Execution Schema, 20 and CI/CD with automatic testing. The pipeline requires a single multi-sample VCF file as input and has a separate workflow for downloading reference datasets and checking their consistency. IBD segments detection workflows of GRAPE can work with both phased and unphased data. As real-world datasets are often heterogeneous and inconsistent, GRAPE incorporates various data preprocessing and quality control (QC) options. GRAPE has a modular architecture that allows switching between tools and adjusting tools' parameters for better control of precision and recall levels. The pipeline also contains a simulation workflow with an in-depth evaluation of pipeline accuracy using simulated and reference data.
GRAPE can work like a standalone version as well as from the dedicated Docker container. We recommend using the containerized version of GRAPE since it has all the dependencies already installed. We published the GRAPE container image in both Docker Hub and Dockstore 21 repositories to satisfy GA4GH 22 standards for sharing Docker-based tools.
As a result, we got a robust, reliable, and easy-to-use tool. Analysis for 10k samples with 600k SNPs requires half an hour, and about 22 hours have been required to process the 100k samples dataset. Precision/recall analysis for GRAPE was performed using simulated datasets. We compared GRAPE with TRIBES, 23 another open-source pipeline for relatedness detection, and showed the advantages of our solution in a sense of the precision/recall metrics.

REVISED Amendments from Version 1
In version 2, we made several changes to address reviewers comments. The main change is the addition of the RaPID IBD inference tool to the GRAPE and its accuracy analysis. We found that RaPID is significantly faster than IBIS when the dataset is phased and its size exceeds 100K. We also tested RaPID-based IBD inference flow on 100K and 500K datasets and provided some computing time estimates for IBIS-based IBD inference flow for 500K datasets. We added Figure 7, which does show that the accuracy of RaPID is similar to accuracy of IBIS except for distinction between 1st and 2nd degree relatives. It requires an additional IBD segment postprocessing to find IBD2 segments and we have not implemented this step yet.
We updated the article and GRAPE source code to ensure that it can process datasets of any size in batches. The memory consumption should not exceed 16GB. Also, we fixed some minor language mistakes and typos.
Any further responses from the reviewers can be found at the end of the article

Methods
This section describes the input data, reference datasets, and the main pipeline steps. The scheme of the pipeline is presented in Figure 1. As input, GRAPE uses a single VCF file containing the genotypes of a set of individuals. Configuration of the pipeline is managed by config.yaml file, or via the parameters of the GRAPE launcher. Reference data should be previously downloaded and stored on a hard drive. After that, the pipeline performs relationship inference accordingly to one of the three possible workflows (described in a corresponding section below). The main pipeline steps are: 1. Downloading of the reference datasets.
2. Quality control and data preprocessing.
3. Application of the relationship inference workflow.
GRAPE also includes simulation workflow to evaluate precision/recall metrics on the simulated data. To simulate artificial pedigrees, Ped-sim tool 24 is used, while unrelated founders for the simulation are taken from the 1000 Genomes Project data. 25 We estimate relationship degrees using pairwise IBD segments count and length. We use IBIS 8 or GERMLINE 7 for the IBD segments inference and ERSA 12 for the estimation of the most likely relationship degree. We included an option to use KING 14 for both IBD segment inference and relationship estimation, however, KING can be used only for close relatives. Also, we added the RaPID 26 option for IBD1 segments inference to compare its accuracy with IBIS and GERMLINE. IBIS is a fast tool that can operate with unphased data. IBIS performs an IBD detection using homozygous single nucleotide polymorphisms (SNPs) only. It breaks the genome into windows of fixed length and searches for homozygous SNP mismatches for each genotype window. If the number of mismatches does not exceed some predefined threshold, the window becomes part of an IBD segment. In contrast, GERMLINE is slower and can work only with phased data, but under some circumstances, it may produce higher precision results. Once the pairwise IBD segments search is completed, relationship degrees are estimated using the ERSA tool. 12 Downloading of the reference datasets GRAPE requires various reference data to perform preprocessing, quality control, phasing, and imputation. Reference data is also required for the simulation workflow. In order to facilitate the collection of reference data, we created a separate workflow that automates this step. It can be run by specifying a reference command to the GRAPE pipeline launcher. The workflow downloads data unpacks it and performs the required post-processing procedures. If phasing and genotype imputation are required, one should also add the additional flags--phase and --impute to the command. It affects the amount of downloaded data. If these flags are specified, the workflow downloads additional reference datasets to make phasing and imputation possible.
There is another option to download all required reference data as a single file. This file is prepared by us and preloaded on our side in the cloud. It can be done by specifying the additional flag --use-bundle to the workflow. This way is faster, since all the post-processing procedures have already been performed. Reference files consist of three main groups.
• Files for the preprocessing. These files include genetic recombination maps for mapping SNPs coordinates from base pairs (bp) to centimorgans (cM); files with the SNPs information from the 1000 Genomes Project for the SNPs quality control; reference genome of hg37 build; liftOver chain file.
• Files for phasing and imputation. Phasing is required as a preliminary step for the GERMLINE tool if input data is unphased. Upon the phasing is done, genotype imputation can be additionally applied. These files are space-demanding and require a considerable amount of post-processing time.
• Files for simulation. These files include phased per-chromosome files from the 1000 Genomes Project; Affymetrix chip data that is used as a source of founders for the simulation; sex-specific recombination maps for better Ped-sim simulation results. 24 Quality control and data preprocessing GRAPE has a versatile and configurable preprocessing workflow. One part of the preprocessing is required and must be performed before the relationship inference workflow. It is launched by the preprocess command of the GRAPE. Along with some necessary technical procedures, preprocessing includes the following steps.
• [Required] SNPs quality control by minor allele frequency (MAF) and the missingness rate. We discovered that blocks of rare SNPs with low MAF values in genotype arrays may produce false positive IBD segments. To address this problem, we filter SNPs by minor allele frequency. We remove SNPs with a MAF value less than 0.02. Additionally, we remove multiallelic SNPs, insertions/deletions, and SNPs with a high missingness rate, because such SNPs are inconsistent with IBD detection tools.
• [Required] Per-sample quality control, using missingness and heterozygosity. Extensive testing revealed that samples with an unusually low level of heterozygosity could produce many false relatives matches among individuals. GRAPE excludes such samples from the analysis and creates a report file with the description of the exclusion reason.
• [Required] Control for strands and SNP IDs mismatches. During this step, GRAPE fixes inconsistencies in strands and reference alleles.
• [Optional] LiftOver from hg38 to hg37. Currently GRAPE uses the hg37 build version of the human genome reference. The pipeline supports input in hg38 and hg37 builds. One should specify the genome build version by the dedicated flag --assembly of the pipeline launcher. If the hg38 build is selected (--assembly hg38), then GRAPE applies the liftOver tool to the input data in order to match the hg37 reference assembly. This parameter should also be specified during the simulation workflow.
• [Optional] Phasing and imputation. GRAPE supports phasing and genotype imputation. GERMLINE IBD detection tool requires phased data. So, if input data is unphased, one should include phasing (--phase flag) into the preprocessing before running the GERMLINE workflow. If input data is highly heterogeneous in a sense of available SNPs positions, we recommend including the imputation procedure as well (--impute).
• [Optional] Removal of imputed SNPs. We found that if the input data is homogeneous in a sense of SNPs positions, the presence of imputed SNPs does not affect the overall IBD detection accuracy of the IBIS tool, but it significantly slows down the overall performance. For this particular case, when input data initially contains a lot of imputed SNPs, we recommend removing them by specifying the --remove-imputation flag to the GRAPE launcher. GRAPE removes all SNPs which are marked with the IMPUTED flag in the input VCF file.

GRAPE workflows
There are three relationship inference workflows implemented in GRAPE. These workflows are activated by the find command of the launcher. Workflow selection is made by the --flow parameter.
1. IBIS + ERSA, --flow ibis. During this workflow, IBD segments detection is performed by IBIS, 8 and estimation of relationship degree is carried out by means of the ERSA algorithm. 12 This is the fastest workflow.
2. IBIS + ERSA & KING, --flow ibis-king. KING 14 is a well-known method for the inference of close relationships. It's fast and can work with unphased data. During this workflow, GRAPE uses the KING tool for the first three degrees of relationships, and the IBIS + ERSA approach for higher-order degrees (see Figure 1). A comparison of evaluation time between IBIS + ERSA and IBIS + ERSA & KING workflows is presented in Figure 3.

GERMLINE + ERSA & KING,
--flow germline-king. The workflow uses GERMLINE for IBD segments detection. KING is used to identify relationships for the first three degrees, and the ERSA algorithm is used for higher order degrees. This workflow was added to GRAPE mainly for the case when input data is already phased and accurately preprocessed.

4.
RaPID + ERSA (BETA) --flow rapid. During this workflow, IBD segments detection is performed by RaPID, 26 and estimation of relationship degree is carried out by means of the ERSA algorithm. It requires input data to be phased, similar to the germline-king workflow. This workflow currently cannot distinguish between parent-offspring and full siblings relationships, because RaPID cannot infer IBD2 segments yet.

Pedigree simulation
We added a simulation workflow into GRAPE to perform a precision/recall analysis of the pipeline. It's accessible by simulating command of the pipeline launcher and incorporates the following steps: (1) pedigree simulation with unrelated founders; here we use the Ped-sim simulation package 24 ; (2) relatedness degrees estimation; (3) comparison between true and estimated degrees. The source dataset for the simulation is taken from CEU (Northern Europeans from Utah) population data of the 1000 Genomes Project. 25 As CEU data consists of trios, we picked no more than one member of each trio as a founder. We also ran GRAPE on selected individuals to remove all cryptic relationships up to the 6th degree. Then, we randomly assigned sex to each individual and used sex-specific genetic maps to take into account the differences in recombination rates between men and women. 27 Results of our precision/recall analysis for the simulated datasets are presented in the corresponding Results section.

IBD segments weighing
Distribution of IBD segments among non-related (background) individuals within a population may be quite heterogeneous. 13 There may exist genome regions with extremely high rates of overall matching, which are not inherited from the recent common ancestors. 28 Instead, these regions are more likely to reflect other demographic factors of the population. The implication is that IBD segments detected in such regions are expected to be less useful for estimating recent relationships. Moreover, such regions are potentially prone to false-positive IBD segments.
GRAPE can use two different approaches to address this issue. The first one is based on the genome regions exclusion mask, wherein some genome regions are completely excluded from the consideration. This approach was proposed by authors of the ERSA algorithm, see Ref. 13. The mask was computed based on whole-genome sequencing data for European individuals. The computed mask is built-in into ERSA 1.0 algorithm and is used by GRAPE by default.
The second approach is based on the so-called IBD segments weighing. This idea reminds one proposed in the Ancestry DNA Matching White Paper 29 (see description of their Timber algorithm there). The key idea is to down-weight the IBD segment, i.e. reduce the IBD segment length, if the segment crosses regions with a high rate of matching. The approach can be briefly described as follows. At first, one should evaluate IBD segments B j È É for some background population. After that, we break chromosomes into windows W i of fixed length and compute total overlap length c i between found IBD segments B j È É and each window W i , c i ¼ P j jW ij j, where W ij is an overlap between the segment B j and the window W i . Obtained overlap lengths are transformed into weights w i which are assigned to the windows, Here f is a weighing function that reflects the following heuristic: if the overlap length c k for a window W k is a relative outlier among all c i f g, then the value of f c k ð Þ is close to 0, otherwise, it's close to 1. Given that weights for all windows are computed, for each IBD segment G one can compute its weighted length G j j w using the formula where jG i j denotes the overlap length between IBD segment G and a window W i . Weighted lengths of IBD segments then are used by the ERSA algorithm (while the ERSA mask is disabled).
GRAPE provides an ability to compute the weight mask from the VCF file with presumably unrelated individuals. It breaks each chromosome into 1cM windows to compute overlaps. After that, GRAPE detects outliers among overlap lengths by means of the minimum covariance determinant (MCD) algorithm, 30 and then determines the outliers upper bound h. This upper bound is used to compute weights, The computed mask can be further used as a parameter for relationship inference workflow (see the --weight-mask parameter). As an example, Figure 2 depicts the weight mask computed for the individuals of East Asian Ancestry taken from the 1000 Genomes Project. 25 To detect IBD segments, IBIS workflow was used with parameters --ibis-seg-len 5, -ibis-min-snp 400. High-matching regions (w i ≈ 0Þ are marked with a blue color. For comparison, ERSA 1.0 masked regions are depicted with the pick hatching.
Operation GRAPE can be run inside a Docker container. This way is recommended. Another option is to run the pipeline from scratch with all the dependencies pre-installed. We successfully ran the pipeline on Ubuntu 18.04 Linux distribution 31 with eight CPUs and 32 GB of RAM to evaluate the performance on huge datasets with up to 100k samples and millions of SNPs. Resource allocation GRAPE can utilize multiple cores. For that, one should specify the number of cores via the --cores parameter. The default number of cores is equal to the total number of available CPUs minus 1.

Execution by scheduler
The pipeline can be run using Funnel, 32 a lightweight task scheduler that implements Task Execution Schema 20 developed by GA4GH. 22 The scheduler can work in various environments, from regular virtual machines to Kubernetes clusters with the support of resource quotas. We provide several examples of the task specifications for Funnel. Each sample represents a JSON file with the task description. These files are available in the GRAPE GitHub repository within the corresponding funnel subfolder.

Performance
To estimate the performance of the GRAPE pipeline we used a machine with eight CPU cores, Intel(R) Xeon(R) CPU E5-2686 v4 @ 2.30 GHz, and 32 GB of RAM. The fastest IBIS + ERSA (--flow ibis) relatedness inference workflow takes about 22 hours to process a 100k individuals dataset, see Figure 3. IBIS tool has a quadratic time complexity with respect to the total number of individuals. The addition of KING (--flow ibis-king) increases the total running time by roughly 50%. Performance analysis confirmed that IBIS is a simple and efficient tool. It allows the pipeline to process hundred of thousands of individuals in a reasonable amount of time. We also evaluated the RaPID + ERSA workflow on the same 100k individuals dataset and its running time was 45 minutes. On 500k dataset, the RaPID + ERSA workflow took 6 hours. However, phasing of datasets of such magnitude typically requires hours or days. The total running time of both workflows including the preprocessing and data preparation stage should be comparable.

Use cases
This section gives examples of the GRAPE pipeline commands which have to be run to infer relationships, or to evaluate precision/recall metrics on a simulated dataset.
Relationship inference with IBIS + ERSA As the first step, reference data must be downloaded. We suppose that reference data is to be stored in the /media/ref directory. Here and below we specify the --real-run flag. Without this flag, GRAPE performs a dry-run.  As the second step, preprocessing is performed. We suppose that the input file is located at/media/input.vcf.gz, and it's in the hg38 build. Input file location specified by the flag --vcf-file. The GRAPE working directory is/media/data. It's specified by the --directory flag.

Listing 2. Preprocessing for the IBIS + ERSA workflow.
The third step is the relationship inference. It is launched with the find command. We use IBIS + ERSA workflow that corresponds to the ibis value of the --flow parameter (default).
Listing 3. Relationship inference with the IBIS + ERSA workflow.
GRAPE has the ability to specify additional parameters to the ERSA and IBIS algorithms to control the sensitivity and the false positive rate. These parameters are described below. Default GRAPE parameters are quite conservative. They provide a low false positive rate and low sensitivity in high (9+) degrees. •

Relationship inference with IBIS + ERSA & KING
The first and second steps are the same as for the previous use case. The third step is launched by specifying the ibis-king flow parameter. For this case, GRAPE performs an estimation of the first three degrees of relationship with KING. Other degrees are estimated with ERSA (see Figure 1). The resulting report of individual relationships for this case contains additional king_degree and kinship columns. KING algorithm has no additional parameters.

Relationship inference with GERMLINE + ERSA & KING
At first, reference data must be downloaded. The GERMLINE tool works with phased data only. So, if input data is unphased, one should download an additional reference dataset to perform phasing. For that purpose, we use a reference panel from the 1000 Genomes Project. This panel takes a considerable amount of disk space ($25 GB) and requires significant time to download. To download this panel one should specify the --phase flag while using the reference command.
The second step is data preprocessing. We use the --phase flag to apply the phasing procedure during this stage.

Evaluation of the IBIS + ERSA workflow on a simulated dataset
To perform simulation, at the first step, one should download the full reference dataset by using the reference command with the --phase and --impute flags enabled. The second step is the simulation workflow. For that, one should use the simulate command of the GRAPE launcher.
Listing 8. Evaluation of the IBIS + ERSA workflow on a simulated dataset.
Along with the parameters for preprocessing, ERSA and IBIS tools, the simulation workflow has two additional parameters listed below.
• --sim-params-file. File with parameters of simulation for the Ped-Sim tool. For more information see. 24 We have prepared several simulation parameters files and stored them in the GRAPE repository on GitHub.
• --sim-samples-file. File with a list of individuals from the 1000 Genomes Project, which are used as founders for the simulations. One can choose ceph_unrelated_all.tsv (unrelated individuals from the CEU population), or all.tsv (all individuals available in 1KGP).
The output files include a list of kinship matches found in the simulated dataset, precision/recall plots, and a confusion matrix to compare the detected degrees of relationship with the true degrees. Detailed information on the computed metrics is presented in the Results section.

Evaluation of GERMLINE + ERSA & KING workflow on a simulated dataset
The first step is to download the reference dataset (see the previous section). The second step is to run the simulate command, specifying germline-king flow and the --phase flag.

Listing 9. Evaluation of the GERMLINE + ERSA & KING workflow on a simulated dataset.
Computation of the IBD segments weighing mask GRAPE has a dedicated command to compute the weight mask and compute-weight-mask. It performs IBD segments detection for the input VCF file and then analyses IBD segments distribution to compute the weight mask. The resulting files consist of a weight mask file in JSON format and a visualization of the mask (see Figure 2).
Listing 10. Computation of the IBD segments weighing mask.
For the example above the resulting files are stored in the/media/background/weight-mask/directory.

Usage of the IBD segments weighing mask
To apply a weight mask during the relatedness detection (find command), one should specify the mask file with the --weight-mask parameter. When used, the ERSA 1.0 exclusion mask is disabled.
Listing 11. Usage of the IBD segments weighing mask.

Results
To test the accuracy and flexibility of the pipeline, we performed extensive testing on both real and simulated datasets. As a sanity check, we took the Allen Ancient DNA Resource (AADR) dataset 33 and made sure that GRAPE does not produce any kinship matches between ancient and present-day individuals.
Next, we have run Khazar origins dataset 34 through the GRAPE. The Khazar dataset contains 1770 samples from 106 Jewish and non-Jewish populations. The dataset contains a significant amount of data from small homogeneous populations. The KING analysis was previously applied to this dataset by the authors 34 to identify close relatives up to the third degree. We applied GRAPE and found 1715 putative relationships. Most of them have a degree of 4+. This result highlights the fact that the total length of IBD segments in small homogeneous populations is several orders of magnitude higher than for the heterogeneous populations in Europe and Asia. 35,36 This is an obvious obstacle to relatedness detection. There is no method known to us to address this issue while using genotypic data obtained from SNP arrays. Whole genome sequencing has the potential to solve this problem since rare mutations should break long IBD segments.
As for the data from genetic testing companies, GRAPE has been successfully applied to the Atlas Biomed database of 100k+ customers, 37 a company that provides direct-to-consumer genetic tests. During this test, GRAPE was proven to be able to handle diverse data obtained from different chips, reference alignments, and other possible data inconsistencies.

Precision/recall analysis on simulated datasets
Finally, GRAPE was evaluated on simulated datasets. We performed the simulation using unrelated founders from 1KGP Affymetrix genotype chip data. This chip contains approximately 900k SNPs. The simulation was carried out with the Ped-Sim package. Using this tool, we produced several pedigree structures with eight generations and a maximum degree of relationships equal to 14. Then we joined all of the pedigrees into one dataset and performed the precision/ recall analysis of the GRAPE pipeline for the different flows. For each degree i of relationships we computed precision and recall metrics: Here TP i ð Þ, FP i ð Þ, FN i ð Þ are the numbers of true positive, false positive, and false negative relationship matches predicted for the degree i. In our analysis, we used non-exact (fuzzy) interval metrics. For the 1st degree, we require an exact match. For the 2nd, 3rd, and 4th degrees, we allow a degree interval of AE1. For example, for the 2nd true degree, we consider a predicted 3rd degree as a true positive match. For the 5th+ degrees, we use the ERSA confidence intervals which are typically 3-4 degrees wide. For 10th+ degrees, these intervals are 6-7 degrees wide. We also plot a confusion matrix for the predicted vs true degrees.
The results of the simulation for the IBIS + ERSA & KING workflow are presented in Figure 4. The following set of parameters was used: --ibis-seg-len 7, --ibis-min-snp 500, --zero-seg-count 0.5, --zero-seg-len 5, --alpha 0.01. The confusion matrix is presented on the right panel in Figure 4. There -1 stands for no-relationship. The pipeline shows recall above 90%+ for degrees from 1 to 5. It detects all relatives with 1-4 degrees. GRAPE found no false positive matches, i.e. it does not find any relationships among unrelated individuals, Recall À1 ð Þ¼1. Precision is above 90% among all detected degrees. We set these parameters as default for this GRAPE workflow. These parameters are quite conservative, i.e. they provide high precision but low sensitivity.
One can relax GRAPE parameters to get more relatedness matches for high degrees. On the other hand, the number of false positive matches increases as well. Figure 5 shows the simulation results for the IBIS + ERSA & KING workflow for slightly relaxed parameters: --ibis-seg-len 5, --ibis-min-snp 400, --zero-seg-count 0.1, --zero-seg-len 5, --alpha 0.01. The number of detected relationships increases significantly for higher (≥8) degrees. But false positive matches arise as well, i.e. Recall À1 ð Þ6 ¼ 1.
Simulation results for the GERMLINE + ERSA & KING workflow are presented in Figure 6. We used the same ERSA parameters as for Figure 4: --zero-seg-count 0.5, --zero-seg-len 5, --alpha 0.01. In comparison to Figure 4, one can see that GERMLINE slightly decreases recall for 6th and 7th degrees, but improves it for higher (≥8) degrees. Precision is above 95% among all detected degrees. No false positive matches were found, i.e. Recall À1 ð Þ¼1. Our experiments showed that, in comparison to IBIS, GERMLINE is better suited for a careful analysis of relatively small cohorts with phased data. The IBIS algorithm is less sensitive, but for segments with lengths above 7 cM, it produces the same results as GERMLINE, while working with unphased data.
We also evaluated the RaPID + ERSA workflow on the same data and found that it has similar accuracy in comparison with IBIS and GERMLINE-based workflows (Figure 7). RaPID works with phased data only, however, it has linear computational complexity with regard to the number of samples. We can recommend RaPID + ERSA workflow for the    distant relationships inference if the data has been already phased. It should work faster and have accuracy comparable with IBIS-based workflows on datasets containing more than 50K samples. Recall and precision are bad for the RaPID workflow because we did not implement a special preprocessing for distinguishing between IBD1 and IBD2 segments. This preprocessing is necessary for the separation of parent-offspring from full siblings relationships. Also, the RaPID + ERSA workflow generates some false positives, particularly for the 6th degree.

Performance of the IBD segments weighing
Our simulation experiments showed that both weighing and exclusion approaches reduce the false-positive rate and significantly improve the overall performance of the pipeline. The weighing mask can be better adapted to specific ancestries, and after additional parameters tuning may slightly outperform ERSA 1.0 approach. In Figure 8, the comparison between the ERSA 1.0 exclusion mask and the GRAPE weight mask is presented for a simulated dataset with the founders of East Asian Ancestry from 1KGP. IBIS workflow was used. The weight mask was taken from Figure 2. Parameter --zero-seg-count was varied while using a weight mask to achieve the same level of precision. One can see, that with approximately the same precision, the weighing approach gives several percentages higher recall.

Comparison with TRIBES
In the end, we compared GRAPE with TRIBES. 23 TRIBES is an earlier open-source pipeline for relatedness detection. The pipeline combines the GERMLINE algorithm for IBD segments detection and the calculation of the genome proportion with zero alleles inferred IBD (IBD0) for each pair to detect the relatedness. If the data is not phased, TRIBES provides the ability to phase data with the EAGLE tool. This part of TRIBES is similar to one of the corresponding GRAPE workflows of IBD segments detection. In contrast to GRAPE, TRIBES estimates degrees of relationship according to expected IBD0 segments proportion ranges. GRAPE uses the ERSA algorithm, which, to our knowledge, is a more advanced approach.
We have run the TRIBES pipeline on the same simulated datasets. Since the simulated datasets contain unphased data, we also applied a built-in phasing procedure from TRIBES. The result of the analysis is presented in Figure 9. TRIBES has demonstrated high detection power for distant relationships up to the 12th degree. Given that many 13+ degree relatives do not share any IBD segments, this is near a theoretical limit. However, TRIBES produces a huge number of false positive matches, see the confusion matrix in the right panel of Figure 9. Since TRIBES lacks options that allow users to control false positive rates by varying pipeline parameters, it becomes a crucial drawback. This obstacle does not allow the TRIBES pipeline to be adapted for applications, where desired precision/recall rate may vary depending on different business or research objectives.

Conclusions
In the current paper, we introduced GRAPE: genomic relatedness detection pipeline. We performed a careful selection of tools and combined various preprocessing steps, IBD segments detection tools, and the relationship estimation algorithm into a single pipeline. One of the possible workflows of the pipeline is based on IBIS tools and can work with unphased data. It's suitable for the analysis of large cohorts in a relatively short time. Using this workflow GRAPE can perform relationship estimation among 100k samples in 22 hours. Another possible workflow is based on the GERMLINE IBD segments detection tool and works with phased data. Our experiments showed that GERMLINE has the most detection power, while the IBIS option is the fastest, easiest to use, and has sufficient accuracy for 1-8 degrees of relationship. Finally, we compared GRAPE with TRIBES, another relatedness detection pipeline. In contrast to GRAPE, TRIBES produce a huge number of false positive matches, requires phased data, and lacks important preprocessing and evaluation options, which makes it impractical for various applications. GRAPE is proved to be a reliable and accurate tool for the analysis of close and distant degrees of kinship. It provides the ability to control a false positive rate, can work with heterogeneous data obtained from various chips, and is ready for production integration.

Data availability
Source data Publicly available datasets were used to test GRAPE. These datasets are available from: • The 1000 Genomes Project 25 ; • URL: https://www.internationalgenome.org/data-portal/data-collection/phase-3.

Open Peer Review
ERSA). The user may want to know if they can run the pipeline in large panels within a reasonable time. The run time of the pipeline using IBD detection (with GERMLINE) was also not included. More recent IBD detection tools (e.g. hapIBD, RaPID and iLASH) can handle very large cohorts efficiently. Incorporating other IBD detection tools might also improve the efficiency of the pipeline regarding the run time and/or memory usage.
Answer: Based on Figure 3 and the fact that IBIS and ERSA algorithms has O(n^2) complexity, when n is the number of samples, we can estimate that running time on n~=500k will be ~20-25 days. This is comparable with time required for phasing a dataset of this size. We are planning to investigate accuracy and performance on hapIBD, RaPID and iLASH for phased datasets, but unique advantage of IBIS-based workflows is that it does not require phasing. We chose GERMLINE as the well-known method with an appropriate license. iLASH does not have an MIT or GPL-like license.

1.
Answer: In the new article version we added beta RaPID support and compared it with IBIS and GERMLINE. RaPID-based workflow took 45 min to run on 100K dataset and 6 hours to run on 500K dataset. But the dataset was already phased.

2.
Robustness against genotyping errors, different marker densities, and different populations have not been studied well. The authors have pointed out the limitation with a small isolated population though.
Answer: IBIS-based workflows detect relatively big IBD segments (more than 5-7 cM). We hope that it will provide some robustness to genotyping error. We also provide a -remove imputation option in the preprocessing workflow, because we found out that some imputed regions produce many false-positive IBD matches.

1.
Answer: IBIS and RaPID-based workflows should not detect many falsepositive IBD segments in the low-density regions, because they constrain both minimum number of SNPs in an IBD match and minimum cM length. However, low-density regions could produce many false-negatives.

2.
3. The design is straightforward, and this pipeline includes different inference strategies that can deal with both phased and unphased data, and it can be run inside a Docker container, which makes it very applicable in real analysis. It is also appreciable that authors devoted lots of times in building this pipeline. I have only three comments: I may understand why to use the King software for close relatedness inference while using IBD based method for distant relatedness inference. It is because large IBD cutoff for example 7cM could lead to loss of IBD segments and bias the kinship estimation. However, when the target population is structured or admixed, the King software may not work as well as IBD based inference.

1.
What should I do if I need to run analysis on a large sample size, for example, UK biobank (n~=500k)? The largest test data in this work is 100k, while some other methods declare they can deal with kinship estimation on much larger sample sizes, see ref 1 .

2.
If I want to find out all potential relatives of one target sample in a reference database, with this pipeline, I need to merge the data first and then detect all relatedness. This is not efficient. An alternative solution is to use query based IBD detection, see ref 2 .