tests/testthat/test_glmm.score.R

In GMMAT: Generalized Linear Mixed Model Association Tests

test_that("cross-sectional id le 400 binomial", {
	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
	pheno <- rbind(example$pheno, example$pheno[1:100, ])
	pheno$id <- 1:500
	pheno$disease[sample(1:500,20)] <- NA
	pheno$age[sample(1:500,20)] <- NA
	pheno$sex[sample(1:500,20)] <- NA
	pheno <- pheno[sample(1:500,450), ]
	pheno <- pheno[pheno$id <= 400, ]
	kins <- example$GRM
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.noselect.1.tmp <- tempfile()
	expect_error(glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1.tmp, ncores = 2), "Error: parallel computing currently not implemented for PLINK binary format genotypes.")
	unlink(obj1.outfile.bed.noselect.1.tmp)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.noselect.1.tmp <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1.tmp, ncores = 2)
	obj1.bgen.noselect.1.tmp <- read.table(obj1.outfile.bgen.noselect.1.tmp, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.1.tmp)
	unlink(obj1.outfile.bgen.noselect.1.tmp)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])

	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.noselect.1.tmp <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1.tmp, ncores = 2)
	obj1.gds.noselect.1.tmp <- read.table(obj1.outfile.gds.noselect.1.tmp, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.1.tmp)
	unlink(obj1.outfile.gds.noselect.1.tmp)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL)), signif(c(0.003804942, 0.986534857)))
	unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt.select.1.tmp <- tempfile()
	expect_error(glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1.tmp, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"), ncores = 2), "Error: parallel computing currently not implemented for plain text format genotypes.")
	unlink(obj1.outfile.txt.select.1.tmp)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)
	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))

	skip_on_cran()

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.003738918, 0.996996766)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("cross-sectional id gt 400 binomial", {
	skip_on_cran()

	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
	pheno <- rbind(example$pheno, example$pheno[1:100, ])
	pheno$id <- 1:500
	pheno$disease[sample(1:500,20)] <- NA
	pheno$age[sample(1:500,20)] <- NA
	pheno$sex[sample(1:500,20)] <- NA
	pheno <- pheno[sample(1:500,450), ]
	kins <- diag(500)
	kins[1:400, 1:400] <- example$GRM
	rownames(kins) <- colnames(kins) <- 1:500
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL)), signif(c(0.00396727, 0.99042091)))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.003967475, 0.992551527)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
	obj1 <- glmmkin(disease ~ age + sex, data = pheno, kins = kins, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

	obj2 <- glmmkin(disease ~ age + sex, data = pheno, kins = NULL, id = "id", family = binomial(link = "logit"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))
	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1, obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("cross-sectional id le 400 gaussian", {
	skip_on_cran()

	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
	pheno <- rbind(example$pheno, example$pheno[1:100, ])
	pheno$id <- 1:500
	pheno$disease[sample(1:500,20)] <- NA
	pheno$age[sample(1:500,20)] <- NA
	pheno$sex[sample(1:500,20)] <- NA
	pheno <- pheno[sample(1:500,450), ]
	pheno <- pheno[pheno$id <= 400, ]
	kins <- example$GRM
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL)), signif(c(0.0316311, 0.9919671)))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.0002293625, 0.9972561847)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))
	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1, obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("cross-sectional id gt 400 gaussian", {
	skip_on_cran()

	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
	pheno <- rbind(example$pheno, example$pheno[1:100, ])
	pheno$id <- 1:500
	pheno$disease[sample(1:500,20)] <- NA
	pheno$age[sample(1:500,20)] <- NA
	pheno$sex[sample(1:500,20)] <- NA
	pheno <- pheno[sample(1:500,450), ]
	kins <- diag(500)
	kins[1:400, 1:400] <- example$GRM
	rownames(kins) <- colnames(kins) <- 1:500
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL)), signif(c(0.0326723, 0.9937123)))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.0003196319, 0.9952539725)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
	obj1 <- glmmkin(trait ~ age + sex, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

	obj2 <- glmmkin(trait ~ age + sex, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"), method = "REML", method.optim = "AI")
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))
	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1, obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("longitudinal repeated measures gaussian", {
	skip_on_cran()

	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
        pheno <- example$pheno2
        kins <- example$GRM
        obj1 <- glmmkin(y.repeated ~ sex, data = pheno, kins = kins, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL)), signif(c(0.01701654, 0.99757462)))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)

        obj2 <- glmmkin(y.repeated ~ sex, data = pheno, kins = NULL, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.01468818, 0.99432063)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
        obj1 <- glmmkin(y.repeated ~ sex, data = pheno, kins = kins, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

        obj2 <- glmmkin(y.repeated ~ sex, data = pheno, kins = NULL, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
        obj1 <- glmmkin(y.repeated ~ sex, data = pheno, kins = kins, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

        obj2 <- glmmkin(y.repeated ~ sex, data = pheno, kins = NULL, id = "id",random.slope = NULL, family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))
	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1, obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("longitudinal random time trend gaussian", {
	skip_on_cran()

	plinkfiles <- strsplit(system.file("extdata", "geno.bed", package = "GMMAT"), ".bed", fixed = TRUE)[[1]]
	bgenfile <- system.file("extdata", "geno.bgen", package = "GMMAT")
	samplefile <- system.file("extdata", "geno.sample", package = "GMMAT")
	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	txtfile <- system.file("extdata", "geno.txt", package = "GMMAT")
	txtfile1 <- system.file("extdata", "geno.txt.gz", package = "GMMAT")
	txtfile2 <- system.file("extdata", "geno.txt.bz2", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(123)
        pheno <- example$pheno2
        kins <- example$GRM
        obj1 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = kins, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.1)
	obj1.bed.noselect.1 <- read.table(obj1.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bed.select.1 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.1)
	obj1.bed.select.1 <- read.table(obj1.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.select.1)
	obj1.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.1)
	obj1.bgen.noselect.1 <- read.table(obj1.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.1)
	obj1.bgen.select.1 <- read.table(obj1.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.select.1)
	expect_equal(obj1.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj1.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(obj1.bed.select.1$PVAL, signif(obj1.gds.select.1$PVAL))
	expect_equal(signif(range(obj1.gds.select.1$PVAL), digits = 5), signif(c(0.02677175, 0.99508480), digits = 5))
	}
	
	obj1.outfile.txt.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.1 <- read.table(obj1.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1$PVAL, obj1.txt.select.1$PVAL)
	obj1.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.1 <- read.table(obj1.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt1.select.1)
	obj1.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.1 <- read.table(obj1.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt2.select.1)

        obj2 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = NULL, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.1)
	obj2.bed.noselect.1 <- read.table(obj2.outfile.bed.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bed.select.1 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.1)
	obj2.bed.select.1 <- read.table(obj2.outfile.bed.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.select.1)
	obj2.outfile.bgen.noselect.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.1)
	obj2.bgen.noselect.1 <- read.table(obj2.outfile.bgen.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.bgen.select.1 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.1)
	obj2.bgen.select.1 <- read.table(obj2.outfile.bgen.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.select.1)
	expect_equal(obj2.bed.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")], obj2.bgen.select.1[, c("SNP","CHR","POS","A1","A2","N","AF","SCORE","VAR","PVAL")])
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(obj2.bed.select.1$PVAL, signif(obj2.gds.select.1$PVAL))
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(0.01366701, 0.98519024)))
	}
	
	obj2.outfile.txt.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.1 <- read.table(obj2.outfile.txt.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1$PVAL, obj2.txt.select.1$PVAL)
	obj2.outfile.txt1.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.1 <- read.table(obj2.outfile.txt1.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt1.select.1)
	obj2.outfile.txt2.select.1 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.1, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.1 <- read.table(obj2.outfile.txt2.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt2.select.1)

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
        obj1 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = kins, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.2)
	obj1.bed.noselect.2 <- read.table(obj1.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.2)
	obj1.outfile.bed.select.2 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.2)
	obj1.bed.select.2 <- read.table(obj1.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.2)
	obj1.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.2)
	obj1.bgen.noselect.2 <- read.table(obj1.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.2)
	obj1.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.2)
	obj1.bgen.select.2 <- read.table(obj1.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)
	}
	
	obj1.outfile.txt.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.2 <- read.table(obj1.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.2)
	obj1.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.2 <- read.table(obj1.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.2)
	obj1.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.2 <- read.table(obj1.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.2)

        obj2 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = NULL, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.2)
	obj2.bed.noselect.2 <- read.table(obj2.outfile.bed.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.2)
	obj2.outfile.bed.select.2 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.2)
	obj2.bed.select.2 <- read.table(obj2.outfile.bed.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.2)
	obj2.outfile.bgen.noselect.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.2)
	obj2.bgen.noselect.2 <- read.table(obj2.outfile.bgen.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.2)
	obj2.outfile.bgen.select.2 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.2)
	obj2.bgen.select.2 <- read.table(obj2.outfile.bgen.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.2)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)
	}
	
	obj2.outfile.txt.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.2 <- read.table(obj2.outfile.txt.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.2)
	obj2.outfile.txt1.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.2 <- read.table(obj2.outfile.txt1.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.2)
	obj2.outfile.txt2.select.2 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.2, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.2 <- read.table(obj2.outfile.txt2.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
        obj1 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = kins, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, outfile = obj1.outfile.bed.noselect.3)
	obj1.bed.noselect.3 <- read.table(obj1.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.noselect.1, obj1.bed.noselect.3)
	obj1.outfile.bed.select.3 <- tempfile()
	glmm.score(obj1, infile = plinkfiles, select = select, outfile = obj1.outfile.bed.select.3)
	obj1.bed.select.3 <- read.table(obj1.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bed.select.1, obj1.bed.select.3)
	obj1.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj1.outfile.bgen.noselect.3)
	obj1.bgen.noselect.3 <- read.table(obj1.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.noselect.1, obj1.bgen.noselect.3)
	obj1.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj1, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj1.outfile.bgen.select.3)
	obj1.bgen.select.3 <- read.table(obj1.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.bgen.select.1, obj1.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)
	}
	
	obj1.outfile.txt.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile, outfile = obj1.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt.select.3 <- read.table(obj1.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt.select.1, obj1.txt.select.3)
	obj1.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile1, outfile = obj1.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt1.select.3 <- read.table(obj1.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt1.select.1, obj1.txt1.select.3)
	obj1.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj1, infile = txtfile2, outfile = obj1.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj1.txt2.select.3 <- read.table(obj1.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.txt2.select.1, obj1.txt2.select.3)

        obj2 <- glmmkin(y.trend ~ sex + time, data = pheno, kins = NULL, id = "id",random.slope = "time", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.bed.noselect.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, outfile = obj2.outfile.bed.noselect.3)
	obj2.bed.noselect.3 <- read.table(obj2.outfile.bed.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.noselect.1, obj2.bed.noselect.3)
	obj2.outfile.bed.select.3 <- tempfile()
	glmm.score(obj2, infile = plinkfiles, select = select, outfile = obj2.outfile.bed.select.3)
	obj2.bed.select.3 <- read.table(obj2.outfile.bed.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bed.select.1, obj2.bed.select.3)
	obj2.outfile.bgen.noselect.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, outfile = obj2.outfile.bgen.noselect.3)
	obj2.bgen.noselect.3 <- read.table(obj2.outfile.bgen.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.noselect.1, obj2.bgen.noselect.3)
	obj2.outfile.bgen.select.3 <- tempfile()
	glmm.score(obj2, infile = bgenfile, BGEN.samplefile = samplefile, select = select, outfile = obj2.outfile.bgen.select.3)
	obj2.bgen.select.3 <- read.table(obj2.outfile.bgen.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.bgen.select.1, obj2.bgen.select.3)
	
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) {
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)
	}
	
	obj2.outfile.txt.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile, outfile = obj2.outfile.txt.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt.select.3 <- read.table(obj2.outfile.txt.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt.select.1, obj2.txt.select.3)
	obj2.outfile.txt1.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile1, outfile = obj2.outfile.txt1.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt1.select.3 <- read.table(obj2.outfile.txt1.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt1.select.1, obj2.txt1.select.3)
	obj2.outfile.txt2.select.3 <- tempfile()
	glmm.score(obj2, infile = txtfile2, outfile = obj2.outfile.txt2.select.3, infile.nrow.skip = 5, infile.ncol.skip = 3, infile.ncol.print = 1:3, select=select, infile.header.print = c("SNP", "Allele1", "Allele2"))
	obj2.txt2.select.3 <- read.table(obj2.outfile.txt2.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.txt2.select.1, obj2.txt2.select.3)

	unlink(c(obj1.outfile.bed.noselect.1, obj1.outfile.bed.select.1, obj1.outfile.bgen.noselect.1, obj1.outfile.bgen.select.1, obj1.outfile.txt.select.1, obj1.outfile.txt1.select.1, obj1.outfile.txt2.select.1))
	unlink(c(obj2.outfile.bed.noselect.1, obj2.outfile.bed.select.1, obj2.outfile.bgen.noselect.1, obj2.outfile.bgen.select.1, obj2.outfile.txt.select.1, obj2.outfile.txt1.select.1, obj2.outfile.txt2.select.1))
	unlink(c(obj1.outfile.bed.noselect.2, obj1.outfile.bed.select.2, obj1.outfile.bgen.noselect.2, obj1.outfile.bgen.select.2, obj1.outfile.txt.select.2, obj1.outfile.txt1.select.2, obj1.outfile.txt2.select.2))
	unlink(c(obj2.outfile.bed.noselect.2, obj2.outfile.bed.select.2, obj2.outfile.bgen.noselect.2, obj2.outfile.bgen.select.2, obj2.outfile.txt.select.2, obj2.outfile.txt1.select.2, obj2.outfile.txt2.select.2))
	unlink(c(obj1.outfile.bed.noselect.3, obj1.outfile.bed.select.3, obj1.outfile.bgen.noselect.3, obj1.outfile.bgen.select.3, obj1.outfile.txt.select.3, obj1.outfile.txt1.select.3, obj1.outfile.txt2.select.3))
	unlink(c(obj2.outfile.bed.noselect.3, obj2.outfile.bed.select.3, obj2.outfile.bgen.noselect.3, obj2.outfile.bgen.select.3, obj2.outfile.txt.select.3, obj2.outfile.txt1.select.3, obj2.outfile.txt2.select.3))
	if(requireNamespace("SeqArray", quietly = TRUE) && requireNamespace("SeqVarTools", quietly = TRUE)) unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1, obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1, obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2, obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2, obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3, obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})

test_that("multiple phenotypes gaussian", {
	skip_on_cran()
	skip_if_not_installed("SeqArray")
	skip_if_not_installed("SeqVarTools")

	gdsfile <- system.file("extdata", "geno.gds", package = "GMMAT")
	data(example)
	suppressWarnings(RNGversion("3.5.0"))
	set.seed(103)
	kins <- example$GRM
	tau1 <- matrix(c(3,0.5,0,0.5,2.5,-0.1,0,-0.1,3),3,3)
	tau2 <- matrix(c(2.5,0.8,0.2,0.8,4.8,-0.1,0.2,-0.1,2.8),3,3)
	kins.chol <- chol(tau1 %x% kins + tau2 %x% diag(400))
	tmp <- as.vector(crossprod(kins.chol, rnorm(1200)))
	x1 <- rnorm(400)
	x2 <- rbinom(400,1,0.5)
	pheno <- data.frame(id = 1:400, x1 = x1, x2 = x2, y1 = 0.5*x1+0.8*x2+tmp[1:400], y2 = x1-0.3*x2+tmp[401:800], y3 = x2+tmp[801:1200])
	obj1 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.1)
	obj1.gds.noselect.1 <- read.table(obj1.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj1.outfile.gds.select.1 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.1)
	obj1.gds.select.1 <- read.table(obj1.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.select.1)
	expect_equal(signif(range(obj1.gds.select.1$PVAL), digits = 5), signif(c(0.009074957, 0.999072499), digits = 5))

	obj2 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.gds.noselect.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.1)
	obj2.gds.noselect.1 <- read.table(obj2.outfile.gds.noselect.1, header = TRUE, as.is = TRUE)
	obj2.outfile.gds.select.1 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.1)
	obj2.gds.select.1 <- read.table(obj2.outfile.gds.select.1, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.select.1)
	expect_equal(signif(range(obj2.gds.select.1$PVAL)), signif(c(2.591209e-05, 9.961137e-01)))

	idx <- sample(nrow(pheno))
	pheno <- pheno[idx, ]
	obj1 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.2)
	obj1.gds.noselect.2 <- read.table(obj1.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.2)
	obj1.outfile.gds.select.2 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.2)
	obj1.gds.select.2 <- read.table(obj1.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.2)

	obj2 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.gds.noselect.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.2)
	obj2.gds.noselect.2 <- read.table(obj2.outfile.gds.noselect.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.2)
	obj2.outfile.gds.select.2 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.2)
	obj2.gds.select.2 <- read.table(obj2.outfile.gds.select.2, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.2)

	idx <- sample(nrow(kins))
	kins <- kins[idx, idx]
	obj1 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = kins, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj1$id_include))
	select[is.na(select)] <- 0
	obj1.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, outfile = obj1.outfile.gds.noselect.3)
	obj1.gds.noselect.3 <- read.table(obj1.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.noselect.1, obj1.gds.noselect.3)
	obj1.outfile.gds.select.3 <- tempfile()
	glmm.score(obj1, infile = gdsfile, select = select, outfile = obj1.outfile.gds.select.3)
	obj1.gds.select.3 <- read.table(obj1.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj1.gds.select.1, obj1.gds.select.3)

	obj2 <- glmmkin(cbind(y1,y2,y3)~x1+x2, data = pheno, kins = NULL, id = "id", family = gaussian(link = "identity"))
	select <- match(1:400, unique(obj2$id_include))
	select[is.na(select)] <- 0
	obj2.outfile.gds.noselect.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, outfile = obj2.outfile.gds.noselect.3)
	obj2.gds.noselect.3 <- read.table(obj2.outfile.gds.noselect.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.noselect.1, obj2.gds.noselect.3)
	obj2.outfile.gds.select.3 <- tempfile()
	glmm.score(obj2, infile = gdsfile, select = select, outfile = obj2.outfile.gds.select.3)
	obj2.gds.select.3 <- read.table(obj2.outfile.gds.select.3, header = TRUE, as.is = TRUE)
	expect_equal(obj2.gds.select.1, obj2.gds.select.3)

	unlink(c(obj1.outfile.gds.noselect.1, obj1.outfile.gds.select.1))
	unlink(c(obj2.outfile.gds.noselect.1, obj2.outfile.gds.select.1))
	unlink(c(obj1.outfile.gds.noselect.2, obj1.outfile.gds.select.2))
	unlink(c(obj2.outfile.gds.noselect.2, obj2.outfile.gds.select.2))
	unlink(c(obj1.outfile.gds.noselect.3, obj1.outfile.gds.select.3))
	unlink(c(obj2.outfile.gds.noselect.3, obj2.outfile.gds.select.3))
})