Serveradmin at 19:25, 23 January 2019

2019-01-23T19:25:07Z

Serveradmin at 19:23, 23 January 2019

2019-01-23T19:23:25Z

Serveradmin: Created page with "NOTITLE ==Mixed models Association Analysis Exercise=="

2019-01-23T17:53:31Z

Created page with "__NOTITLE__ ==Mixed models Association Analysis Exercise=="

New page

__NOTITLE__

==Mixed models Association Analysis Exercise==

← Older revision		Revision as of 19:25, 23 January 2019
Line 89:		Line 89:
	hist(blupsgenos$BLUP[blupsgenos$combo==21], main="genotype=21", xlab="BLUP")		hist(blupsgenos$BLUP[blupsgenos$combo==21], main="genotype=21", xlab="BLUP")
	hist(blupsgenos$BLUP[blupsgenos$combo==22], main="genotype=22", xlab="BLUP")		hist(blupsgenos$BLUP[blupsgenos$combo==22], main="genotype=22", xlab="BLUP")
		+	q()
		+	gcta64 --bfile quantfamdata --blup-snp GCTAblups.indi.blp --out GCTAsnpblups
		+	R
		+	snpblups<-read.table("GCTAsnpblups.snp.blp", header=F)
		+	names(snpblups)=c("SNP","Allele","snpBLUP","snpResidual")
		+	head(snpblups)
		+	res2<-read.table("FLMMresults", header=T)
		+	new2<-data.frame(res2$SNP, res2$Chromosome, res2$Position, res2$Pvalue)
		+	names(new2)<-c("SNP", "CHR", "BP", "P")
		+	head(new2)
		+	mergeblup=merge(snpblups,new2, by="SNP",sort=F)
		+	head(mergeblup)
		+	par(mfrow=c(1,1))
		+	plot(mergeblup$snpBLUP,-log10(mergeblup$P))

← Older revision		Revision as of 19:23, 23 January 2019
Line 2:		Line 2:

	==Mixed models Association Analysis Exercise==		==Mixed models Association Analysis Exercise==
		+
		+	head FLMMresults
		+	more topsnps.txt
		+	plink --bfile quantfamdata --extract topsnps.txt --recodeA --out topsnpsFLMM
		+	R
		+	fam<-read.table("quantfamdata.fam", header=F)
		+	pheno=data.frame(fam[,1:2],fam[,6])
		+	write.table(pheno,file="phenos.txt",col.names=F,row.names=F,quote=F)
		+	gcta64 --mlma --bfile quantfamdata --pheno phenos.txt --out GCTAresults
		+	R
		+	source("qqmanHJCupdated.R")
		+	res3<-read.table("GCTAresults.mlma", header=T)
		+	head(res3)
		+	chi<-(qchisq(1-res3$p,1))
		+	lambda=median(chi)/0.456
		+	lambda
		+	new3<-data.frame(res3$SNP, res3$Chr, res3$bp, res3$p)
		+	names(new3)<-c("SNP", "CHR", "BP", "P")
		+	head(new3)
		+	qq(new3$P)
		+	manhattan(new3, pch=20, suggestiveline=F, genomewideline=F, ymin=2, cex.x.axis=0.65, colors=c("black","dodgerblue"), cex=0.5)
		+	res2<-read.table("FLMMresults", header=T)
		+	new2<-data.frame(res2$SNP, res2$Chromosome, res2$Position, res2$Pvalue)
		+	names(new2)<-c("SNP", "CHR", "BP", "P")
		+	merged=merge(new3,new2, by="SNP", sort=F)
		+	head(res2)
		+	head(new2)
		+	head(new3)
		+	head(merged)
		+	plot(-log10(merged$P.x),-log10(merged$P.y))
		+	abline(0,1)
		+	q()
		+	gcta64 --bfile quantfamdata --autosome --make-grm-bin --out GCTAgrm
		+	gcta64 --reml --grm-bin GCTAgrm --pheno phenos.txt --out GCTAherit
		+	gcta64 --bfile quantfamdata --chr 1 --make-grm-bin --out GCTAgrmchr1
		+	gcta64 --reml --grm-bin GCTAgrmchr1 --pheno phenos.txt --out GCTAheritchr1
		+	gcta64 --bfile quantfamdata --chr 2 --make-grm-bin --out GCTAgrmchr2
		+	gcta64 --reml --grm-bin GCTAgrmchr2 --pheno phenos.txt --out GCTAheritchr2
		+	gcta64 --bfile quantfamdata --chr 6 --make-grm-bin --out GCTAgrmchr6
		+	gcta64 --reml --grm-bin GCTAgrmchr6 --pheno phenos.txt --out GCTAheritchr6
		+	gcta64 --bfile quantfamdata --chr 12 --make-grm-bin --out GCTAgrmchr12
		+	gcta64 --reml --grm-bin GCTAgrmchr12 --pheno phenos.txt --out GCTAheritchr12
		+	gcta64 --bfile quantfamdata --chr 1 --make-grm-bin --out GCTAgrmchr1
		+	gcta64 --bfile quantfamdata --chr 2 --make-grm-bin --out GCTAgrmchr2
		+	gcta64 --bfile quantfamdata --chr 3 --make-grm-bin --out GCTAgrmchr3
		+	gcta64 --bfile quantfamdata --chr 4 --make-grm-bin --out GCTAgrmchr4
		+	gcta64 --bfile quantfamdata --chr 5 --make-grm-bin --out GCTAgrmchr5
		+	gcta64 --bfile quantfamdata --chr 6 --make-grm-bin --out GCTAgrmchr6
		+	gcta64 --bfile quantfamdata --chr 7 --make-grm-bin --out GCTAgrmchr7
		+	gcta64 --bfile quantfamdata --chr 8 --make-grm-bin --out GCTAgrmchr8
		+	gcta64 --bfile quantfamdata --chr 9 --make-grm-bin --out GCTAgrmchr9
		+	gcta64 --bfile quantfamdata --chr 10 --make-grm-bin --out GCTAgrmchr10
		+	gcta64 --bfile quantfamdata --chr 11 --make-grm-bin --out GCTAgrmchr11
		+	gcta64 --bfile quantfamdata --chr 12 --make-grm-bin --out GCTAgrmchr12
		+	gcta64 --bfile quantfamdata --chr 13 --make-grm-bin --out GCTAgrmchr13
		+	gcta64 --bfile quantfamdata --chr 14 --make-grm-bin --out GCTAgrmchr14
		+	gcta64 --bfile quantfamdata --chr 15 --make-grm-bin --out GCTAgrmchr15
		+	gcta64 --bfile quantfamdata --chr 16 --make-grm-bin --out GCTAgrmchr16
		+	gcta64 --bfile quantfamdata --chr 17 --make-grm-bin --out GCTAgrmchr17
		+	gcta64 --bfile quantfamdata --chr 18 --make-grm-bin --out GCTAgrmchr18
		+	gcta64 --bfile quantfamdata --chr 19 --make-grm-bin --out GCTAgrmchr19
		+	gcta64 --bfile quantfamdata --chr 20 --make-grm-bin --out GCTAgrmchr20
		+	gcta64 --bfile quantfamdata --chr 21 --make-grm-bin --out GCTAgrmchr21
		+	gcta64 --bfile quantfamdata --chr 22 --make-grm-bin --out GCTAgrmchr22
		+	gcta64 --reml --mgrm-bin multipleGRMs.txt --pheno phenos.txt --out GCTAherit22GRMs
		+	gcta64 --mgrm-bin multipleGRMsnot6.txt --make-grm --out GCTAgrmnot6
		+	gcta64 --reml --mgrm-bin multipleGRMs6andnot6.txt --pheno phenos.txt --out GCTAherit6andnot6
		+	gcta64 --reml --reml-pred-rand --grm-bin GCTAgrm --pheno phenos.txt --out GCTAblups
		+	R
		+	blups<-read.table("GCTAblups.indi.blp", header=F)
		+	names(blups)=c("FID","IID","intV1","BLUP", "intV2", "Residual")
		+	head(blups)
		+	hist(blups$BLUP)
		+	topgenos<-read.table("topsnpsFLMM.raw",header=T)
		+	blupsgenos=data.frame(blups, topgenos)
		+	blupsgenos$combo=10*blupsgenos$rs9271209_G+blupsgenos$rs233722_G
		+	head(blupsgenos)
		+	par(mfrow=c(3,3))
		+	hist(blupsgenos$BLUP[blupsgenos$combo==0], main="genotype=0", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==1], main="genotype=1", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==2], main="genotype=2", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==10], main="genotype=10", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==11], main="genotype=11", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==12], main="genotype=12", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==20], main="genotype=20", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==21], main="genotype=21", xlab="BLUP")
		+	hist(blupsgenos$BLUP[blupsgenos$combo==22], main="genotype=22", xlab="BLUP")

Mixed models Association Analysis Exercise - Revision history

Serveradmin at 19:25, 23 January 2019

Serveradmin at 19:23, 23 January 2019

Serveradmin: Created page with "__NOTITLE__ ==Mixed models Association Analysis Exercise=="

Serveradmin: Created page with "NOTITLE ==Mixed models Association Analysis Exercise=="