Asked by a user via email:

I'm interested in using GENESIS package to infer kinship for an ad-mixed population. When I follow steps on as on http://bioconductor.org/packages/release/bioc/vignettes/GENESIS/inst/doc/pcair.R I get error for subscript out of bounds converting king output to matrix.

Could you help me understand what's wrong?

library("SNPRelate",lib.loc="~/R_LIB")
library("GENESIS",lib.loc="~/R_LIB")
library("GWASTools",lib.loc="~/R_LIB")

bed<-"selected.bed"
fam<-"selected.fam" #10 people in sample data
bim<-"selected.bim"

snpgdsBED2GDS(bed, fam, bim, "9040.gds")
snpgdsSummary("9040.gds")

# read in GDS data
geno <- GdsGenotypeReader(filename = "9040.gds")
# create a GenotypeData class object
genoData <- GenotypeData(geno)

iids <- getScanID(genoData)
head(iids)

#Kinship using KING
# king -b selected.bed  --kinship --prefix  small_data
#
KINGmat <- king2mat(file.kin0 = "small_data.kin0", iids = iids, file.kin="small_data.kin")

##############
##############
#Error: 
#KINGmat <- king2mat(file.kin0 = "small_data.kin0", iids = iids, file.kin="small_data.kin")
#Reading .kin0 file...
#Reading .kin file...
#Determining Unique Individual IDs from KING Output...
#Checking Provided Individual IDs
#Adding Kinship Entries from .kin0 file...
#Adding Kinship Entries from .kin file...
#Error in `[<-`(`*tmp*`, id1[r], id2[r], value = tmp1$Kinship[r]) :
#  subscript out of bounds

sessionInfo()
R version 3.4.2 (2017-09-28)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Debian GNU/Linux 9 (stretch)

Matrix products: default
BLAS: /mnt/mfs/cluster/bin/R-3.4/lib/libRblas.so
LAPACK: /mnt/mfs/cluster/bin/R-3.4/lib/libRlapack.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
 [9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods
[8] base

other attached packages:
[1] GWASTools_1.24.1    Biobase_2.38.0      BiocGenerics_0.24.0
[4] GENESIS_2.10.0      SNPRelate_1.14.0    gdsfmt_1.14.1

loaded via a namespace (and not attached):
 [1] zoo_1.8-1              purrr_0.2.4            DNAcopy_1.52.0
 [4] splines_3.4.2          lattice_0.20-35        GWASExactHW_1.01
 [7] stats4_3.4.2           mgcv_1.8-23            blob_1.1.1
[10] survival_2.42-3        rlang_0.2.0            pillar_1.2.2
[13] glue_1.2.0             DBI_0.8                bit64_0.9-7
[16] bindrcpp_0.2.2         GenomeInfoDbData_1.0.0 bindr_0.1.1
[19] zlibbioc_1.24.0        MatrixModels_0.4-1     Biostrings_2.46.0
[22] memoise_1.1.0          SeqArray_1.18.2        IRanges_2.12.0
[25] SparseM_1.77           lmtest_0.9-36          GenomeInfoDb_1.14.0
[28] quantreg_5.35          Rcpp_0.12.16           quantsmooth_1.44.0
[31] S4Vectors_0.16.0       graph_1.56.0           XVector_0.18.0
[34] bit_1.1-12             digest_0.6.15          dplyr_0.7.4
[37] GenomicRanges_1.30.3   grid_3.4.2             bitops_1.0-6
[40] magrittr_1.5           sandwich_2.4-0         tibble_1.4.2
[43] RCurl_1.95-4.10        RSQLite_2.1.0          mice_2.46.0
[46] tidyr_0.8.0            pkgconfig_2.0.1        MASS_7.3-50
[49] Matrix_1.2-14          assertthat_0.2.0       logistf_1.22
[52] R6_2.2.2               SeqVarTools_1.16.1     rpart_4.1-13
[55] nnet_7.3-12            nlme_3.1-137           compiler_3.4.2

I think something is wrong with the format of your input files from KING. You can compare them with the example files in the expected format provided with the GENESIS package: "MXL_ASW.kin0" and "MXL_ASW.kin" in the "extdata" directory.

Alternatively, you can use the SNPRelate package to run KING in R:

library(SNPRelate)
gdsfile <- system.file("extdata", "HapMap_ASW_MXL_geno.gds", package="GENESIS")
gds <- openfn.gds(gdsfile)
king <- snpgdsIBDKING(gds)
KINGmat <- king$kinship
dimnames(KINGmat) <- list(king$sample.id, king$sample.id)