#This script belongs to the paper of Ketelaar, Portelli, Dijk et al. #Entitled 'Phenotypic and functional translation of IL33 genetics in asthma', submitted to JACI
#Below example scripts performing eQTL analyses using R statistics #Core Team R. R: A language and environment for statistical computing. 2013; Available at: http://www.R-project.org/.
Due to ethical restrictions we can only share code and summary statistics, no raw input data files
set working directory
setwd("H:/x/x/x")
Load snpdata
source('./snpdata.load.R')
Load and combine covariate and phenotype data:
data<- read.table("covariatedata.txt",sep="\t", header=TRUE) pheno<- read.table("pheno.txt",sep="\t",quote="", header=TRUE) ##save(data,pheno,file="data_pheno.Rdata") ##load ("data_pheno.Rdata") ## fast loading data eqtl.index<- which(pheno$IncludeForGWAS==TRUE) ## select samples data<-data[,eqtl.index]; pheno<-pheno[eqtl.index,]
Get individual covariates
smoking<- pheno$Smoking_status smoking[which(is.na(smoking)==TRUE)]<-"3" ## put missing value as "3" gender <- pheno$Gender age <- as.numeric(pheno$Age) diseases <- pheno$Diseases tid<- pheno$TID cov<- as.data.frame(cbind(smoking,gender,age,diseases,tid)) temp<- as.character(pheno$Sample_ID) centers<- apply(t(temp),2, function (x) strsplit(x,"_")1[1]) rm(temp)
Perform regression analyses:
lm.regression<- function(data,snp,cov,pc) {
gender <- cov$gender smoking <- cov$smoking age <- as.numeric(cov$age) diseases <- cov$diseases snp<- unlist(snp) lm.reg <- lm(data~snp+gender+smoking+age+diseases+pc) out<- summary(lm.reg) beta<- out$coefficients[2,"Estimate"] sd <- out$coefficients[2,"Std. Error"] pvalue <- out$coefficients[2, "Pr(>|t|)"] output<- c(beta,sd,pvalue) return(output) }