Process use case Daniel: use R code from within PLINK

[richelbilderbeek]

From Daniel:

It's a built-in feature:https://www.cog-genomics.org/plink/1.9/rserve
(It's plink 1 exclusive, though)
I've been having a lot of trouble with it, though. Plink didn't feel like cooperating.

library(VGAM)
library(foreach)

fit_tobit <- function(pheno, geno, covariates, lower) {
  model <- vglm(pheno ~ geno + covariates, tobit(Lower=lower))
  result <- c(length(pheno), coef(summary(model))[3,]) # Genotype is first predictor and there are 2 intercepts
  c(length(result), result)
}

Rplink <- function(PHENO, GENO, CLUSTER, COVAR) {
  foreach(i=1:ncol(GENO), .combine="c") %do% fit_tobit(PHENO, GENO[,i], COVAR, 175)
}