2016-july-berlin-commands - Revision history

Serveradmin: /* Annotation exercise */

2016-07-08T12:27:27Z

‎Annotation exercise

Serveradmin: /* Annotation exercise */

2016-07-08T12:25:07Z

‎Annotation exercise

Serveradmin: /* Annotation exercise */

2016-07-08T12:19:06Z

‎Annotation exercise

Serveradmin: /* VAT exercise */

2016-06-15T19:14:56Z

‎VAT exercise

Serveradmin: /* RV-TDT exercise */

2016-06-15T00:10:35Z

‎RV-TDT exercise

Serveradmin at 23:59, 14 June 2016

2016-06-14T23:59:47Z

Serveradmin: /* SEQPower exercise */

2016-06-13T20:24:35Z

‎SEQPower exercise

Serveradmin: /* VAT exercise */

2016-06-13T14:52:04Z

‎VAT exercise

Serveradmin: /* VAT exercise */

2016-06-13T14:51:21Z

‎VAT exercise

Serveradmin: Created page with "NOTITLE ==Annotation exercise== mkdir APOC3 && cd APOC3 vtools init APOC3 vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19 vtools select..."

2016-06-13T14:34:05Z

Created page with "__NOTITLE__ ==Annotation exercise== mkdir APOC3 && cd APOC3 vtools init APOC3 vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19 vtools select..."

New page

__NOTITLE__

==Annotation exercise==

mkdir APOC3 && cd APOC3
vtools init APOC3
vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19
vtools select variant --count
vtools use refGene
vtools show fields
vtools select variant --output chr pos ref alt refGene.name2
vtools use ../data/dbNSFP.DB
vtools show fields
vtools select variant 'id="rs76353203" or id="rs140621530" or id="rs147210663" or id="rs121918381" or id="rs121918382"' -t 5var
vtools select 5var --output chr pos ref alt id dbNSFP.GERP_RS dbNSFP.phyloP100way_vertebrate dbNSFP.phastCons100way_vertebrate
vtools select 5var --output chr pos ref alt dbNSFP.FATHMM_pred dbNSFP.FATHMM_score
vtools select 5var --output chr pos ref alt dbNSFP.Polyphen2_HDIV_pred dbNSFP.Polyphen2_HDIV_score dbNSFP.Polyphen2_HVAR_pred dbNSFP.Polyphen2_HVAR_score
vtools select 5var --output chr pos ref alt dbNSFP.LRT_pred dbNSFP.LRT_score
vtools select 5var --output chr pos ref alt dbNSFP.MutationTaster_pred
vtools select 5var --output chr pos ref alt dbNSFP.PROVEAN_pred dbNSFP.PROVEAN_score
vtools select 5var --output chr pos ref alt dbNSFP.SIFT_pred dbNSFP.SIFT_score
vtools select 5var --output chr pos ref alt dbNSFP.CADD_raw dbNSFP.CADD_phred
==GATK exercise==

java -jar GenomeAnalysisTK.jar --help java -Xmx200m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T UnifiedGenotyper -glm BOTH -I child.bam -I father.bam -I mother.bam -L trio.intervals -o ug.vcf
less -S ug.vcf java -Xmx200m -jar GenomeAnalysisTK.jar -T VariantEval -R human_g1k_v37.fa -L trio.intervals -D dbsnp.vcf -eval ug.vcf -o ug.vcf.eval less -S ug.vcf.eval java -Xmx200m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T SelectVariants -V ug.vcf -mv -mvq 0 -o ug.dnm.vcf -ped trio.ped java -Xmx500m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T HaplotypeCaller -I child.bam -I father.bam -I mother.bam -minPruning 4 -L 1:199325670-199325672 -mergeVariantsViaLD -o hc.vcf
==PLINKSEQ exercise==
Data analysis:

pseq help pseq help all pseq myproj new-project --resources hg19 pseq myproj load-vcf --vcf CEU.exon.2010_03.genotypes.hg19.vcf.gz YRI.exon.2010_03.genotypes.hg19.vcf.gz pseq myproj load-pheno --file phenotype.phe pseq myproj v-view | head pseq myproj i-view | head pseq myproj summary pseq myproj var-summary pseq myproj ind-summary pseq myproj loc-summary pseq myproj ref-summary pseq myproj seq-summary pseq myproj file-summary pseq myproj meta-summary pseq myproj v-stats pseq myproj i-stats | head pseq myproj tag-file --id 1 --name CEU pseq myproj tag-file --id 2 --name YRI pseq myproj var-summary pseq myproj v-freq | head pseq myproj v-freq --mask file=CEU | head pseq myproj v-freq --mask file=YRI | head pseq myproj v-view --mask any.filter.ex | head pseq myproj v-view --mask any.filter.ex | wc -l pseq myproj v-view --mask any.filter | wc -l pseq myproj var-set --group pass --mask any.filter.ex pseq myproj var-summary pseq myproj var-set --group pass_DP15 --mask include="DP>14" var=pass pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10 --mask geno=DP:ge:11 var=pass_DP15 pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_CEU --mask file=CEU var=pass_DP15_DPgeno10 pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_CEU_HWE --mask hwe=5.7e-7:1 var=pass_DP15_DPgeno10_CEU pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_CEU_HWE_MAFgt05 --mask maf=0.05:0.5 var=pass_DP15_DPgeno10_CEU_HWE pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_CEU_HWE_MAFlt01 --mask "mac=1 maf=0.01" var=pass_DP15_DPgeno10_CEU_HWE pseq myproj var-summary pseq myproj glm --phenotype BMI --covar SEX --mask var=pass_DP15_DPgeno10_CEU_HWE_MAFgt05 > SNV_CEU.result head SNV_CEU.result cat SNV_CEU.result | awk '{if(FNR==1) print $0; if(NR>1) print $0 | "sort -k9"}' | grep -v "NA\s\+NA\s\+NA" | head pseq myproj assoc --tests skat --phenotype BMI --covar SEX --mask var=pass_DP15_DPgeno10_CEU_HWE_MAFlt01 loc.group=refseq > SKAT_CEU.result pseq myproj assoc --tests skat --phenotype BMI --covar SEX --mask include="DP>14" geno=DP:ge:11 file=CEU hwe=5.7e-7:1 "mac=1 maf=0.01" loc.group=refseq > SKAT_CEU.result head -20 SKAT_CEU.result cat SKAT_CEU.result | grep SKAT | grep -v "P=NA" | sort -k6 | head -15

Exercise analyzing YRI samples:

pseq myproj var-set --group pass_DP15_DPgeno10_YRI --mask file=YRI var=pass_DP15_DPgeno10 pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_YRI_HWE --mask hwe=5.7e-7:1 var=pass_DP15_DPgeno10_YRI pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_YRI_HWE_MAFgt05 --mask maf=0.05:0.5 var=pass_DP15_DPgeno10_YRI_HWE pseq myproj var-summary pseq myproj var-set --group pass_DP15_DPgeno10_YRI_HWE_MAFlt01 --mask "mac=1 maf=0.01" var=pass_DP15_DPgeno10_YRI_HWE pseq myproj var-summary pseq myproj glm --phenotype BMI --covar SEX --mask var=pass_DP15_DPgeno10_YRI_HWE_MAFgt05 > SNV_YRI.result head SNV_YRI.result cat SNV_YRI.result | awk '{if(FNR==1) print $0; if(NR>1) print $0 | "sort -k9"}' | grep -v "NA\s\+NA\s\+NA" | head pseq myproj assoc --tests skat --phenotype BMI --covar SEX --mask include="DP>14" geno=DP:ge:11 file=YRI hwe=5.7e-7:1 "mac=1 maf=0.01" loc.group=refseq > SKAT_YRI.result head -20 SKAT_YRI.result cat SKAT_YRI.result | grep SKAT | grep -v "P=NA" | sort -k6 | head -15

==Popgen exercise R programs==
popgen_drift.q

############################################################# ### Michael Nothnagel, michael.nothnagel@uni-koeln.de ### ### Simulating genetic drift ### #############################################################
 N.s = c(10, 100, 1000, 10000) n.gen = 50 n.rep = 100
# === Simulate genetic drift === # for (n in N.s) { freqs = matrix(NA, ncol=n.gen+1, nrow=n.rep) for (i in 1:n.rep) { alleles = c(rep(0,n/2), rep(1,n/2)) freqs[i,1] = sum(alleles==1) / length(alleles) for (j in 1:n.gen) { alleles = sample (alleles, length(alleles), replace=T) freqs[i,j+1] = sum(alleles==1) / length(alleles) } } assign (paste("freqs.n", n, sep=""), freqs) }
summary(freqs.n10 [,n.gen+1]) summary(freqs.n100 [,n.gen+1]) summary(freqs.n1000 [,n.gen+1]) summary(freqs.n10000[,n.gen+1])
# === Graph allele frequency changes === # pdf("drift_plot.pdf", paper="special", height=4*2, width=4*2, onefile=F) split.screen(c(2,2)) screen(1) plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency", main="N=10") lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA") freqs = get(paste("freqs.n", 10 , sep="")) for (i in 1:n.rep) { lines(0:n.gen,freqs[i,], col="#44AAAA") } screen(2) plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency", main="N=100") lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA") freqs = get(paste("freqs.n", 100 , sep="")) for (i in 1:n.rep) { lines(0:n.gen,freqs[i,], col="#44AAAA") } screen(3) plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency", main="N=1000") lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA") freqs = get(paste("freqs.n", 1000 , sep="")) for (i in 1:n.rep) { lines(0:n.gen,freqs[i,], col="#44AAAA") } screen(4) plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency", main="N=10000") lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA") freqs = get(paste("freqs.n", 10000, sep="")) for (i in 1:n.rep) { lines(0:n.gen,freqs[i,], col="#44AAAA") } close.screen(all.screens=T) dev.off()
# === Graph allele frequency after 50 generations === # pdf("drift_hist.pdf", paper="special", height=4*2, width=4*2, onefile=F) split.screen(c(2,2)) screen(1) hist(freqs.n10 [,n.gen+1], main="N=10", xlab="Allele frequency", xlim=c(0,1)) screen(2) hist(freqs.n100 [,n.gen+1], main="N=100", xlab="Allele frequency", xlim=c(0,1)) screen(3) hist(freqs.n1000 [,n.gen+1], main="N=1000", xlab="Allele frequency", xlim=c(0,1)) screen(4) hist(freqs.n10000[,n.gen+1], main="N=10000", xlab="Allele frequency", xlim=c(0,1)) close.screen(all.screens=T) dev.off()

popgen.selection.q

############################################################# ### Michael Nothnagel, michael.nothnagel@uni-koeln.de ### ### Calculation allele frequency changes due to selection ### #############################################################
 s.s = c(0.001, 0.01, 0.1) h = 0.5 n.gen = 100
# === Calculate allele frequency changes === # for (s in s.s) { freqs = rep(NA, n.gen+1) af = 0.5 freqs[1] = 1-af for (j in 1:n.gen) { omega = 1 - 2*af*(1-af)*h*s - (1-af)*(1-af)*s f.het = (1-h*s)*2*af*(1-af)/omega f.hom = af*af/omega af = f.hom + f.het/2 freqs[j+1] = 1-af } assign (paste("freqs.s", s, sep=""), freqs) }
# === Report allele frequencies after 100 generations === # for (s in s.s) { cat("s=") ; cat(s) ; cat(": ") freqs = get(paste("freqs.s", s, sep="")) cat(freqs[n.gen+1]) ; cat("\n") }
# === Graph allele frequency changes === # pdf("selection_plot.pdf", paper="special", height=4*2, width=4*2, onefile=F) plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency") lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA") for (s in s.s) { freqs = get(paste("freqs.s", s, sep="")) lines(0:n.gen,freqs, col="#44AAAA") } dev.off()
== Regression exercise==
In R:

load("dbp.R") ls() dbp[1:5,] # result.snp12 = glm (affection ~ rs1112, family=binomial("logit"), data=dbp) print (result.snp12) print ( class (result.snp12) ) print ( summary(result.snp12) ) # dev.geno = anova (result.snp12, test="Chi") lrt.pvalue = pchisq(dev.geno[dim(dev.geno)[1],"Deviance"], df=2, ncp=0, FALSE) print ( lrt.pvalue ) # print ( summary(result.snp12)$coefficients ) snp.beta = summary(result.snp12)$coefficients[2:3,1] print ( snp.beta ) print ( exp(snp.beta) ) ci = confint (result.snp12) print (ci) print ( exp(ci) ) # snp.data = dbp[,c("affection", "rs1112")] summary(snp.data) snp.data[,"rs1112"] summary(snp.data) # result.all = glm (affection ~ rs1112, family=binomial("logit"), data=snp.data) dev.all = anova (result.all, test="Chi") summary(result.all) print(dev.all) # snp.data = dbp[,c("affection", "trait","sex", "age", "rs1112", "rs1117")] summary(snp.data) snp.data[,"rs1112"] snp.data[,"rs1117"]
# result.adj = glm (affection ~ sex + rs1112 , family=binomial("logit"), data=snp.data) summary(result.adj) # result.adj = glm (affection ~ age + rs1112 , family=binomial("logit"), data=snp.data) summary(result.adj) # result.adj = glm (affection ~ sex + age + rs1112, family=binomial("logit"), data=snp.data) summary(result.adj) # result.adj = glm (affection ~ rs1117 + rs1112, family=binomial("logit"), data=snp.data) summary(result.adj) anova (result.adj, test="Chi") result.adj = glm (affection ~ rs1112 + rs1117, family=binomial("logit"), data=snp.data) summary(result.adj) anova (result.adj, test="Chi") # result.adj = lm (trait ~ rs1112, data=snp.data) summary(result.adj) result.adj = lm (trait ~ sex + rs1112, data=snp.data) summary(result.adj) # result.inter = glm (affection ~ sex * rs1112, family=binomial("logit"), data=snp.data) summary(result.inter) result.inter = glm (affection ~ age * rs1112, family=binomial("logit"), data=snp.data) summary(result.inter) # result.inter = glm (affection ~ rs1112 * rs1117, family=binomial("logit"), data=snp.data) summary(result.inter) # q()
==SEQPower exercise==
spower -h spower LOGIT -h spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --method "CFisher --name CMC" -r 100 -j 4 -l 1 -o exercise
 spower show exercise.csv spower show exercise.csv power* spower show exercise.loci.csv spower show exercise.loci.csv maf spower show tests spower show test SKAT
 spower LOGIT Kryukov2009European1800.sfs --def_rare 0.01 --def_neutral -0.00001 0.00001 --moi A --proportion_detrimental 1 --proportion_protective 0 --OR_rare_detrimental 1.5 --OR_common_detrimental 1 --baseline_effect 0.01 --sample_size 1000 --p1 0.5 --limit 1 --alpha 0.05 --method "KBAC --name K1 --mafupper 0.01 --maflower 0 --alternative 1 --moi additive --permutations 1000 --adaptive 0.1" --replicates 1000 --jobs 4 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.2 --ORmax_rare_detrimental 3.0 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.2 --ORmax_rare_detrimental 3.0 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower show exercise.loci.csv effect*
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.2 --ORmax_rare_detrimental 3.0 --proportion_detrimental 0.8 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --missing_sites 0.05 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --missing_sites 0.05 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --missing_low_maf 0.000125 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --missing_low_maf 0.000125 --method CFisher -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --method "CFisher --alternative 1 --name CMC" "KBAC --permutations 1000 --alternative 1" "WSSRankTest --alternative 1 --name WSS" "VTtest --alternative 1 --permutations 1000" "SKAT disease" -r 100 -j 4 -l 1 -o exercise spower LNR Kryukov2009European1800.sfs --sample_size 1000 --meanshift_rare_detrimental 0.2 --method "CollapseQt --name CMC --alternative 2" -r 100 -j 4 -l 1 -o exercise
 spower LNR Kryukov2009European1800.sfs --sample_size 1000 --meanshift_rare_detrimental 0.2 --meanshiftmax_rare_detrimental 0.5 --method "CollapseQt --alternative 2" -r 100 -j 4 -l 1 -o exercise
 spower ELNR Kryukov2009European1800.sfs --sample_size 1000 --meanshift_rare_detrimental 0.2 --QT_thresholds 0.4 0.6 --method "CollapseQt --alternative 2" -r 100 -j 4 -l 1 -o exercise
 spower ELNR Kryukov2009European1800.sfs --sample_size 1000 --p1 0.5 --meanshift_rare_detrimental 0.5 --QT_thresholds 0.4 0.6 --method "CollapseQt --alternative 2" -r 100 -j 4 -l 1 -o exercise
 spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --method "GroupWrite ExerciseSimulation" -j 4 -o exercise -v1
 spower show exercise.SEQPowerDB LOGIT method power title --condition "where power between 0.25 and 0.95" 
 for i in 1 1.5 2 2.5 3 3.5 4; do spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental $i --method "CFisher --name CMC$i" --title FixedOR$i -r 100 -j 4 -l 1 -o exercise2 done
==VAT exercise==
vtools -h
vtools init VATDemo
vtools import *.vcf.gz --var_info DP filter --geno_info DP_geno --build hg18 -j1
vtools liftover hg19
head phenotypes.csv
vtools phenotype --from_file phenotypes.csv --delimiter ","
vtools show project
vtools show tables
vtools show table variant
vtools show samples
vtools show genotypes
vtools show fields
vtools select variant --count
vtools show genotypes > GenotypeSummary.txt
head GenotypeSummary.txt
vtools output variant "max(DP)" "min(DP)" "avg(DP)" "stdev(DP)" "lower_quartile(DP)" "upper_quartile(DP)" --header
vtools select variant "filter='PASS'" --count
vtools select variant "filter='PASS'" -o "max(DP)" "min(DP)" "avg(DP)" "stdev(DP)" "lower_quartile(DP)" "upper_quartile(DP)" --header
vtools update variant --from_stat 'total=#(GT)' 'num=#(alt)' 'het=#(het)' 'hom=#(hom)' 'other=#(other)' 'minDP=min(DP_geno)' 'maxDP=max(DP_geno)' 'meanDP=avg(DP_geno)' 'maf=maf()'
vtools show fields
vtools show table variant
vtools update variant --from_stat 'totalGD10=#(GT)' 'numGD10=#(alt)' 'hetGD10=#(het)' 'homGD10=#(hom)' 'otherGD10=#(other)' 'mafGD10=maf()' --genotypes "DP_geno > 10"
vtools show fields
vtools show table variant
vtools output variant chr pos maf mafGD10 --header --limit 20
vtools phenotype --set "RACE=0" --samples "filename like 'YRI%'"
vtools phenotype --set "RACE=1" --samples "filename like 'CEU%'"
vtools show samples --limit 10
vtools update variant --from_stat 'CEU_mafGD10=maf()' --genotypes 'DP_geno>10' --samples "RACE=1"
vtools update variant --from_stat 'YRI_mafGD10=maf()' --genotypes 'DP_geno>10' --samples "RACE=0"
vtools output variant chr pos mafGD10 CEU_mafGD10 YRI_mafGD10 --header --limit 10
vtools phenotype --from_stat 'CEU_totalGD10=#(GT)' 'CEU_numGD10=#(alt)' --genotypes 'DP_geno>10' --samples "RACE=1"
vtools phenotype --from_stat 'YRI_totalGD10=#(GT)' 'YRI_numGD10=#(alt)' --genotypes 'DP_geno>10' --samples "RACE=0"
vtools phenotype --output sample_name CEU_totalGD10 CEU_numGD10 YRI_totalGD10 YRI_numGD10 --header
vtools select variant 'maf>=0.01' -t variant_MAFge01 'Variants that have MAF >= 0.01'
vtools show tables
vtools execute KING --var_table variant_MAFge01
vtools_report plot_pheno_fields KING_MDS1 KING_MDS2 RACE --dot KING.mds.race.pdf --discrete_color Dark2
vtools_report plot_pheno_fields KING_MDS1 KING_MDS2 panel --dot KING.mds.panel.pdf --discrete_color Dark2
vtools execute ANNOVAR geneanno
vtools output variant chr pos ref alt mut_type --limit 20 --header
vtools_report trans_ratio variant -n num
vtools_report trans_ratio variant -n numGD10
vtools select variant "DP<15" -t to_remove
vtools show tables
vtools remove variants to_remove -v0
vtools show tables
vtools remove genotypes "DP_geno<10" -v0 vtools select variant "mut_type like 'non%' or mut_type like 'stop%' or region_type='splicing'" -t v_funct vtools show tables vtools show samples --limit 5 vtools select variant --samples "RACE=1" -t CEU mkdir -p ceu
cd ceu
vtools init ceu --parent ../ --variants CEU --samples "RACE=1" --build hg19
vtools show project
vtools select variant "CEU_mafGD10>=0.05" -t common_ceu
vtools select v_funct "CEU_mafGD10<0.01" -t rare_ceu vtools use refGene vtools show annotation refGene vtools associate -h vtools show tests vtools show test LinRegBurden
vtools associate common_ceu BMI --covariate SEX -m "LinRegBurden --alternative 2" -j1 --to_db EA_CV > EA_CV.asso.res
grep -i error *.log
less EA_CV.asso.res
sort -g -k7 EA_CV.asso.res | head
vtools show fields
vtools associate rare_ceu BMI --covariate SEX -m "LinRegBurden --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV.asso.res
grep -i error *.log | tail -22
less EA_RV.asso.res
sort -g -k6 EA_RV.asso.res | head
vtools associate rare_ceu BMI --covariate SEX -m "VariableThresholdsQt --alternative 2 -p 100000 --adaptive 0.0005" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_VT.asso.res
grep -i error *.log | tail -22
less EA_RV_VT.asso.res
sort -g -k6 EA_RV_VT.asso.res | head
vtools select rare_ceu "refGene.name2='ABCC1'" -o chr pos ref alt CEU_mafGD10 numGD10 mut_type --header
vtools_report plot_association qq -o QQRV -b --label_top 2 -f 6 < EA_RV.asso.res
vtools_report plot_association manhattan -o MHRV -b --label_top 5 --color Dark2 --chrom_prefix None -f 6 < EA_RV.asso.res vtools associate rare_ceu BMI --covariate SEX KING_MDS1 KING_MDS2 -m "LinRegBurden --name RVMDS2 --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_MDS2.asso.res
vtools_report plot_association qq -o QQRV_MDS2 -b --label_top 2 -f 6 < EA_RV_MDS2.asso.res cd .. vtools select variant --samples "RACE=0" -t YRI mkdir -p yri cd yri
vtools init yri --parent ../ --variants YRI --samples "RACE=0" --build hg19
vtools select variant "YRI_mafGD10>=0.05" -t common_yri
vtools select v_funct "YRI_mafGD10<0.01" -t rare_yri vtools use refGene vtools associate common_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -j1 --to_db YA_CV > YA_CV.asso.res
vtools associate rare_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -g refGene.name2 -j1 --to_db YA_RV > YA_RV.asso.res
vtools associate rare_yri BMI --covariate SEX -m "VariableThresholdsQt --alternative 2 -p 100000 --adaptive 0.0005" -g refGene.name2 -j1 --to_db YA_RV > YA_RV_VT.asso.res
cd ..
vtools_report meta_analysis ceu/EA_RV_VT.asso.res yri/YA_RV_VT.asso.res --beta 5 --pval 6 --se 7 -n 2 --link 1 > META_RV_VT.asso.res
cut -f1,3 META_RV_VT.asso.res | head

@@ Line 56: / Line 56: @@
   load("dbp.R")<br data-attributes="%20/" />ls()<br data-attributes="%20/" />dbp[1:5,]<br data-attributes="%20/" />#<br data-attributes="%20/" />result.snp12 = glm (affection ~ rs1112, family=binomial("logit"), data=dbp)<br data-attributes="%20/" />print (result.snp12)<br data-attributes="%20/" />print ( class (result.snp12) )<br data-attributes="%20/" />print ( summary(result.snp12) )<br data-attributes="%20/" />#<br data-attributes="%20/" />dev.geno = anova (result.snp12, test="Chi")<br data-attributes="%20/" />lrt.pvalue = pchisq(dev.geno[dim(dev.geno)[1],"Deviance"],<br data-attributes="%20/" /> df=2, ncp=0, FALSE)<br data-attributes="%20/" />print ( lrt.pvalue )<br data-attributes="%20/" />#<br data-attributes="%20/" />print ( summary(result.snp12)$coefficients )<br data-attributes="%20/" />snp.beta = summary(result.snp12)$coefficients[2:3,1]<br data-attributes="%20/" />print ( snp.beta )<br data-attributes="%20/" />print ( exp(snp.beta) )<br data-attributes="%20/" />ci = confint (result.snp12)<br data-attributes="%20/" />print (ci)<br data-attributes="%20/" />print ( exp(ci) )<br data-attributes="%20/" />#<br data-attributes="%20/" />snp.data = dbp[,c("affection", "rs1112")]<br data-attributes="%20/" />summary(snp.data)<br data-attributes="%20/" />snp.data[,"rs1112"]<br />summary(snp.data)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.all = glm (affection ~ rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />dev.all = anova (result.all, test="Chi")<br data-attributes="%20/" />summary(result.all)<br data-attributes="%20/" />print(dev.all)<br data-attributes="%20/" />#<br data-attributes="%20/" />snp.data = dbp[,c("affection", "trait","sex", "age", "rs1112", "rs1117")]<br data-attributes="%20/" />summary(snp.data)<br data-attributes="%20/" />snp.data[,"rs1112"]<br />snp.data[,"rs1117"]
   #<br data-attributes="%20/" />result.adj = glm (affection ~ sex + rs1112 , family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ age + rs1112 , family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ sex + age + rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ rs1117 + rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />anova (result.adj, test="Chi")<br data-attributes="%20/" />result.adj = glm (affection ~ rs1112 + rs1117, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />anova (result.adj, test="Chi")<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = lm (trait ~ rs1112, data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />result.adj = lm (trait ~ sex + rs1112, data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.inter = glm (affection ~ sex * rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />result.inter = glm (affection ~ age * rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.inter = glm (affection ~ rs1112 * rs1117, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />#<br data-attributes="%20/" />q()
 ==SEQPower exercise==
   spower -h<br data-attributes="%20/" />spower LOGIT -h<br data-attributes="%20/" />spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --method "CFisher --name CMC" -r 100 -j 4 -l 1 -o exercise

@@ Line 21: / Line 21: @@
   vtools select 5var --output chr pos ref alt dbNSFP.SIFT_pred dbNSFP.SIFT_score
   vtools select 5var --output chr pos ref alt dbNSFP.CADD_raw dbNSFP.CADD_phred
-==GATK exercise==
+==GATK/IGV exercise==
   java -jar GenomeAnalysisTK.jar --help<br data-attributes="%20/" />java -Xmx200m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T UnifiedGenotyper -glm BOTH -I child.bam -I father.bam -I mother.bam -L trio.intervals -o ug.vcf
   less -S ug.vcf<br data-attributes="%20/" />java -Xmx200m -jar GenomeAnalysisTK.jar -T VariantEval -R human_g1k_v37.fa -L trio.intervals -D dbsnp.vcf -eval ug.vcf -o ug.vcf.eval<br data-attributes="%20/" />less -S ug.vcf.eval<br data-attributes="%20/" />java -Xmx200m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T SelectVariants -V ug.vcf -mv -mvq 0 -o ug.dnm.vcf -ped trio.ped<br data-attributes="%20/" />java -Xmx500m -jar GenomeAnalysisTK.jar -R human_g1k_v37.fa -T HaplotypeCaller -I child.bam -I father.bam -I mother.bam -minPruning 4 -L 1:199325670-199325672 -mergeVariantsViaLD -o hc.vcf
-==PLINKSEQ exercise==
+==PSEQ exercise==
 Data analysis:
@@ Line 51: / Line 51: @@
   # === Report allele frequencies after 100 generations === #<br data-attributes="%20/" />for (s in s.s) {<br data-attributes="%20/" /> cat("s=") ; cat(s) ; cat(": ")<br data-attributes="%20/" /> freqs = get(paste("freqs.s", s, sep=""))<br data-attributes="%20/" /> cat(freqs[n.gen+1]) ; cat("\n")<br data-attributes="%20/" /> }
   # === Graph allele frequency changes === #<br data-attributes="%20/" />pdf("selection_plot.pdf", paper="special", height=4*2, width=4*2, onefile=F)<br data-attributes="%20/" /> plot(x=0, y=0, type="n", xlim=c(0,n.gen), ylim=c(0,1), xlab="Generation", ylab="Allele frequency")<br data-attributes="%20/" /> lines (c(0,n.gen), rep(0.5, 2), lty=3, col="#AAAAAA")<br data-attributes="%20/" /> for (s in s.s) {<br data-attributes="%20/" /> freqs = get(paste("freqs.s", s, sep=""))<br data-attributes="%20/" /> lines(0:n.gen,freqs, col="#44AAAA")<br data-attributes="%20/" /> }<br data-attributes="%20/" />dev.off()
-==<br data-attributes="%20class=%22Apple-interchange-newline%22">Regression exercise==
+==Regression exercise==
 In R:

@@ Line 56: / Line 56: @@
   load("dbp.R")<br data-attributes="%20/" />ls()<br data-attributes="%20/" />dbp[1:5,]<br data-attributes="%20/" />#<br data-attributes="%20/" />result.snp12 = glm (affection ~ rs1112, family=binomial("logit"), data=dbp)<br data-attributes="%20/" />print (result.snp12)<br data-attributes="%20/" />print ( class (result.snp12) )<br data-attributes="%20/" />print ( summary(result.snp12) )<br data-attributes="%20/" />#<br data-attributes="%20/" />dev.geno = anova (result.snp12, test="Chi")<br data-attributes="%20/" />lrt.pvalue = pchisq(dev.geno[dim(dev.geno)[1],"Deviance"],<br data-attributes="%20/" /> df=2, ncp=0, FALSE)<br data-attributes="%20/" />print ( lrt.pvalue )<br data-attributes="%20/" />#<br data-attributes="%20/" />print ( summary(result.snp12)$coefficients )<br data-attributes="%20/" />snp.beta = summary(result.snp12)$coefficients[2:3,1]<br data-attributes="%20/" />print ( snp.beta )<br data-attributes="%20/" />print ( exp(snp.beta) )<br data-attributes="%20/" />ci = confint (result.snp12)<br data-attributes="%20/" />print (ci)<br data-attributes="%20/" />print ( exp(ci) )<br data-attributes="%20/" />#<br data-attributes="%20/" />snp.data = dbp[,c("affection", "rs1112")]<br data-attributes="%20/" />summary(snp.data)<br data-attributes="%20/" />snp.data[,"rs1112"]<br />summary(snp.data)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.all = glm (affection ~ rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />dev.all = anova (result.all, test="Chi")<br data-attributes="%20/" />summary(result.all)<br data-attributes="%20/" />print(dev.all)<br data-attributes="%20/" />#<br data-attributes="%20/" />snp.data = dbp[,c("affection", "trait","sex", "age", "rs1112", "rs1117")]<br data-attributes="%20/" />summary(snp.data)<br data-attributes="%20/" />snp.data[,"rs1112"]<br />snp.data[,"rs1117"]
   #<br data-attributes="%20/" />result.adj = glm (affection ~ sex + rs1112 , family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ age + rs1112 , family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ sex + age + rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = glm (affection ~ rs1117 + rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />anova (result.adj, test="Chi")<br data-attributes="%20/" />result.adj = glm (affection ~ rs1112 + rs1117, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />anova (result.adj, test="Chi")<br data-attributes="%20/" />#<br data-attributes="%20/" />result.adj = lm (trait ~ rs1112, data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />result.adj = lm (trait ~ sex + rs1112, data=snp.data)<br data-attributes="%20/" />summary(result.adj)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.inter = glm (affection ~ sex * rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />result.inter = glm (affection ~ age * rs1112, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />#<br data-attributes="%20/" />result.inter = glm (affection ~ rs1112 * rs1117, family=binomial("logit"),<br data-attributes="%20/" /> data=snp.data)<br data-attributes="%20/" />summary(result.inter)<br data-attributes="%20/" />#<br data-attributes="%20/" />q()
 ==SEQPower exercise==
   spower -h<br data-attributes="%20/" />spower LOGIT -h<br data-attributes="%20/" />spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental 1.5 --method "CFisher --name CMC" -r 100 -j 4 -l 1 -o exercise
@@ Line 75: / Line 89: @@
   <br data-attributes="%20/" />spower show exercise.SEQPowerDB LOGIT method power title --condition "where power between 0.25 and 0.95"<br data-attributes="%20/" />
   <br data-attributes="%20/" />for i in 1 1.5 2 2.5 3 3.5 4; do<br data-attributes="%20/" />spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental $i --method "CFisher --name CMC$i" --title FixedOR$i -r 100 -j 4 -l 1 -o exercise2<br data-attributes="%20/" />done
 ==VAT exercise==

@@ Line 166: / Line 166: @@
   vtools_report meta_analysis ceu/EA_RV_VT.asso.res yri/YA_RV_VT.asso.res --beta 5 --pval 6 --se 7 -n 2 --link 1 > META_RV_VT.asso.res
   cut -f1,3 META_RV_VT.asso.res | head==VMT==
   ##################################
   # Example 1: Autosomal recessive #

@@ Line 2: / Line 2: @@
 ==Annotation exercise==
   vtools init APOC3
   vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19

@@ Line 4: / Line 4: @@
   mkdir APOC3 &amp;&amp; cd APOC3
   vtools init APOC3
-  vtools import ../APOC3.vcf --format ../vcf.fmt --build hg19
+  vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19
   vtools select variant --count
   vtools use refGene
   vtools show fields
   vtools select variant --output chr pos ref alt refGene.name2
-  vtools use ../dbNSFP.DB
+  vtools use ../data/dbNSFP.DB
   vtools show fields
   vtools select variant 'id="rs76353203" or id="rs140621530" or id="rs147210663" or id="rs121918381" or id="rs121918382"' -t 5var

@@ Line 77: / Line 77: @@
 ==VAT exercise==
   vtools -h
   vtools init VATDemo
@@ Line 125: / Line 124: @@
   vtools remove variants to_remove -v0
   vtools show tables
-  vtools remove genotypes "DP_geno<10" -v0 vtools select variant "mut_type like 'non%' or mut_type like 'stop%' or region_type='splicing'" -t v_funct
+  vtools remove genotypes "DP_geno<10" -v0
   vtools show tables
   vtools show samples --limit 5
@@ Line 134: / Line 134: @@
   vtools show project
   vtools select variant "CEU_mafGD10>=0.05" -t common_ceu
-  vtools select v_funct "CEU_mafGD10<0.01" -t rare_ceu vtools use refGene
+  vtools select v_funct "CEU_mafGD10<0.01" -t rare_ceu
   vtools show annotation refGene
   vtools associate -h
@@ Line 154: / Line 155: @@
   vtools select rare_ceu "refGene.name2='ABCC1'" -o chr pos ref alt CEU_mafGD10 numGD10 mut_type --header
   vtools_report plot_association qq -o QQRV -b --label_top 2 -f 6 < EA_RV.asso.res
-  vtools_report plot_association manhattan -o MHRV -b --label_top 5 --color Dark2 --chrom_prefix None -f 6 < EA_RV.asso.res vtools associate rare_ceu BMI --covariate SEX KING_MDS1 KING_MDS2 -m "LinRegBurden --name RVMDS2 --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_MDS2.asso.res
+  vtools_report plot_association manhattan -o MHRV -b --label_top 5 --color Dark2 --chrom_prefix None -f 6 < EA_RV.asso.res
-  vtools_report plot_association qq -o QQRV_MDS2 -b --label_top 2 -f 6 < EA_RV_MDS2.asso.res cd ..
+ vtools associate rare_ceu BMI --covariate SEX KING_MDS1 KING_MDS2 -m "LinRegBurden --name RVMDS2 --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_MDS2.asso.res
   vtools select variant --samples "RACE=0" -t YRI
   mkdir -p yri; cd yri
   vtools init yri --parent ../ --variants YRI --samples "RACE=0" --build hg19
   vtools select variant "YRI_mafGD10>=0.05" -t common_yri
-  vtools select v_funct "YRI_mafGD10<0.01" -t rare_yri vtools use refGene
+  vtools select v_funct "YRI_mafGD10<0.01" -t rare_yri
   vtools associate common_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -j1 --to_db YA_CV > YA_CV.asso.res
   vtools associate rare_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -g refGene.name2 -j1 --to_db YA_RV > YA_RV.asso.res
@@ Line 167: / Line 171: @@
   vtools_report meta_analysis ceu/EA_RV_VT.asso.res yri/YA_RV_VT.asso.res --beta 5 --pval 6 --se 7 -n 2 --link 1 > META_RV_VT.asso.res
   cut -f1,3 META_RV_VT.asso.res | head==VMT==

@@ Line 76: / Line 76: @@
   <br data-attributes="%20/" />for i in 1 1.5 2 2.5 3 3.5 4; do<br data-attributes="%20/" />spower LOGIT Kryukov2009European1800.sfs --sample_size 1000 --OR_rare_detrimental $i --method "CFisher --name CMC$i" --title FixedOR$i -r 100 -j 4 -l 1 -o exercise2<br data-attributes="%20/" />done
 ==VAT exercise==
   vtools -h
   vtools init VATDemo
@@ Line 123: / Line 124: @@
   vtools remove variants to_remove -v0
   vtools show tables
-  vtools remove genotypes "DP_geno<10" -v0 <br />vtools select variant "mut_type like 'non%' or mut_type like 'stop%' or region_type='splicing'" -t v_funct <br />vtools show tables <br />vtools show samples --limit 5 <br />vtools select variant --samples "RACE=1" -t CEU <br />mkdir -p ceu
+  vtools remove genotypes "DP_geno<10" -v0 vtools select variant "mut_type like 'non%' or mut_type like 'stop%' or region_type='splicing'" -t v_funct
   cd ceu
   vtools init ceu --parent ../ --variants CEU --samples "RACE=1" --build hg19
   vtools show project
   vtools select variant "CEU_mafGD10>=0.05" -t common_ceu
-  vtools select v_funct "CEU_mafGD10<0.01" -t rare_ceu <br />vtools use refGene <br />vtools show annotation refGene <br />vtools associate -h <br />vtools show tests <br />vtools show test LinRegBurden
+  vtools select v_funct "CEU_mafGD10<0.01" -t rare_ceu vtools use refGene
   vtools associate common_ceu BMI --covariate SEX -m "LinRegBurden --alternative 2" -j1 --to_db EA_CV > EA_CV.asso.res
   grep -i error *.log
@@ Line 135: / Line 144: @@
   vtools show fields
   vtools associate rare_ceu BMI --covariate SEX -m "LinRegBurden --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV.asso.res
-  grep -i error *.log | tail -22
+  grep -i error *.log | tail -10
   less EA_RV.asso.res
   sort -g -k6 EA_RV.asso.res | head
   vtools associate rare_ceu BMI --covariate SEX -m "VariableThresholdsQt --alternative 2 -p 100000 --adaptive 0.0005" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_VT.asso.res
-  grep -i error *.log | tail -22
+  grep -i error *.log | tail -10
   less EA_RV_VT.asso.res
   sort -g -k6 EA_RV_VT.asso.res | head
   vtools select rare_ceu "refGene.name2='ABCC1'" -o chr pos ref alt CEU_mafGD10 numGD10 mut_type --header
   vtools_report plot_association qq -o QQRV -b --label_top 2 -f 6 < EA_RV.asso.res
-  vtools_report plot_association manhattan -o MHRV -b --label_top 5 --color Dark2 --chrom_prefix None -f 6 < EA_RV.asso.res <br />vtools associate rare_ceu BMI --covariate SEX KING_MDS1 KING_MDS2 -m "LinRegBurden --name RVMDS2 --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_MDS2.asso.res
+  vtools_report plot_association manhattan -o MHRV -b --label_top 5 --color Dark2 --chrom_prefix None -f 6 < EA_RV.asso.res vtools associate rare_ceu BMI --covariate SEX KING_MDS1 KING_MDS2 -m "LinRegBurden --name RVMDS2 --alternative 2" -g refGene.name2 -j1 --to_db EA_RV > EA_RV_MDS2.asso.res
-  vtools_report plot_association qq -o QQRV_MDS2 -b --label_top 2 -f 6 < EA_RV_MDS2.asso.res <br />cd .. <br />vtools select variant --samples "RACE=0" -t YRI<br />mkdir -p yri<br />cd yri
+  vtools_report plot_association qq -o QQRV_MDS2 -b --label_top 2 -f 6 < EA_RV_MDS2.asso.res cd ..
   vtools init yri --parent ../ --variants YRI --samples "RACE=0" --build hg19
   vtools select variant "YRI_mafGD10>=0.05" -t common_yri
-  vtools select v_funct "YRI_mafGD10<0.01" -t rare_yri <br />vtools use refGene  vtools associate common_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -j1 --to_db YA_CV > YA_CV.asso.res
+  vtools select v_funct "YRI_mafGD10<0.01" -t rare_yri vtools use refGene
   vtools associate rare_yri BMI --covariate SEX -m "LinRegBurden --alternative 2" -g refGene.name2 -j1 --to_db YA_RV > YA_RV.asso.res
   vtools associate rare_yri BMI --covariate SEX -m "VariableThresholdsQt --alternative 2 -p 100000 --adaptive 0.0005" -g refGene.name2 -j1 --to_db YA_RV > YA_RV_VT.asso.res
   cd ..
   vtools_report meta_analysis ceu/EA_RV_VT.asso.res yri/YA_RV_VT.asso.res --beta 5 --pval 6 --se 7 -n 2 --link 1 > META_RV_VT.asso.res
-  cut -f1,3 META_RV_VT.asso.res | head
+  cut -f1,3 META_RV_VT.asso.res | head==VMT==

2016-july-berlin-commands - Revision history

Serveradmin: /* Annotation exercise */

Serveradmin: /* Annotation exercise */

Serveradmin: /* Annotation exercise */

Serveradmin: /* VAT exercise */

Serveradmin: /* RV-TDT exercise */

Serveradmin at 23:59, 14 June 2016

Serveradmin: /* SEQPower exercise */

Serveradmin: /* VAT exercise */

Serveradmin: /* VAT exercise */

Serveradmin: Created page with "__NOTITLE__ ==Annotation exercise== mkdir APOC3 && cd APOC3 vtools init APOC3 vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19 vtools select..."

Serveradmin: Created page with "NOTITLE ==Annotation exercise== mkdir APOC3 && cd APOC3 vtools init APOC3 vtools import ../data/APOC3.vcf --format ../data/vcf.fmt --build hg19 vtools select..."