require(devtools)
install_bitbucket("kindlychung/txtutils")
# install_bitbucket("kindlychung/manqq2")
# install_github("kindlychung/gmailR")
# install_bitbucket("kindlychung/autoGmail")
install_bitbucket("kindlychung/collr2")
library(collr2)
setwd("~/Desktop/mmp3")
plinkArgs = getPlinkParam(allow_no_sex = "", pheno = "mmp13.phe", pheno_name = "Page", covar = "mmp13.phe", covar_name = "Sex,Cage", linear = "hide-covar")
initGwasArgs = getPlinkParam(maf=.01, hwe=1e-4, allow_no_sex = "", pheno = "mmp13.phe", pheno_name = "Page", assoc= "")
crt1("~/Desktop/mmp3", "test4", plinkArgs, initGwasArgs, FALSE, 0.1, 5)
debug(permute)
pv = permute(plinkArgs = plinkArgs, initGwasArgs = initGwasArgs, phenoFileOrig = "mmp13.phe", pFilter = .2, nShiftMax = 10, n=10)
sort(pv)
require(collr2)
setwd("/media/data1/kaiyin/RS123_1KG/")
plinkArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "EasilyBurn",
covar = "RS123.1kg.pheno/dermatology.csv",
covar_name = "sex,age,SC,SunProtect",
logistic = "hide-covar",
one = ""
)
initGwasArgs = getPlinkParam(
allow_no_sex = "",
missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "EasilyBurn",
one = "",
assoc=""
)
crt1(
wDir = "/media/data1/kaiyin/RS123_1KG",
taskName = "easilyburn2_s200_p0.01",
plinkArgs = plinkArgs,
initGwasArgs = initGwasArgs,
initGwas = TRUE,
pFilter = 1e-2,
nMaxShift = 200
)
require(collr2)
setwd("/media/data1/kaiyin/RS123_1KG")
plinkArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "AK",
covar = "RS123.1kg.pheno/dermatology.csv",
covar_name = "sex,age,SC",
linear = "hide-covar"
)
initGwasArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "AK",
assoc = ""
)
crt1(
wDir = "/media/data1/kaiyin/RS123_1KG",
taskName = "AK_s200_p0.01",
plinkArgs = plinkArgs,
initGwasArgs = initGwasArgs,
initGwas = TRUE,
pFilter = 1e-2,
nMaxShift = 200
)
#############################
## sagging
#############################
require(collr2)
setwd("/media/data1/NL_AU_UK/NL_AU_UK_sagging_replicate/chr5frag/")
plinkArgs = getPlinkParam(
allow_no_sex = "",
linear = "hide-covar"
)
initGwasArgs = getPlinkParam(
allow_no_sex = "",
assoc=""
)
routine2(
pheno = "NL_AU_UK_sagging1.pheno",
pheno_name = "sag",
covar_name = "sex,age",
plinkArgs = plinkArgs,
initGwasArgs = initGwasArgs,
pFilter = 1e-1,
nMaxShift = 400
)
require(collr2)
setwd("/media/data1/kaiyin/RS123_1KG")
plinkArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "Sagging",
covar = "RS123.1kg.pheno/dermatology.csv",
covar_name = "sex,age,SC",
linear = "hide-covar"
)
initGwasArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "Sagging",
assoc = ""
)
crt1(
wDir = "/media/data1/kaiyin/RS123_1KG",
taskName = "Sagging_s200_p0.01",
plinkArgs = plinkArgs,
initGwasArgs = initGwasArgs,
initGwas = TRUE,
pFilter = 1e-2,
nMaxShift = 200
)
require(collr2)
setwd("/media/data1/kaiyin/RS123_1KG")
plinkArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "Log_PS",
covar = "RS123.1kg.pheno/dermatology.csv",
covar_name = "sex,age,SC",
linear = "hide-covar"
)
initGwasArgs = getPlinkParam(
allow_no_sex = "", missing_phenotype = 9999,
pheno = "RS123.1kg.pheno/dermatology.csv",
pheno_name = "Log_PS",
assoc = ""
)
crt1(
wDir = "/media/data1/kaiyin/RS123_1KG",
taskName = "Log_PS_s200_p0.01",
plinkArgs = plinkArgs,
initGwasArgs = initGwasArgs,
initGwas = TRUE,
pFilter = 1e-2,
nMaxShift = 200
)
require(collr2)
setwd("/media/data1/kaiyin/RS123_1KG/")
globalMinVec = numeric(0)
phenoFileOrig = "RS123.1kg.pheno/dermatology.csv"
phenoFile = "~/Desktop/shuffle_pheno.csv"
for(i in 1:100) {
i = 1
taskName = sprintf("easilyburn2_permutation_%d", i)
shuffleLines(phenoFileOrig, phenoFile, "EasilyBurn")
plinkArgs = getPlinkParam(allow_no_sex = "", missing_phenotype = 9999, pheno = phenoFile, covar = phenoFile, covar_name = "sex,age,SC,SunProtect", logistic = "hide-covar", one = "", pheno_name = "EasilyBurn")
initGwasArgs = getPlinkParam(allow_no_sex = "", missing_phenotype = 9999, pheno = phenoFile, one = "", pheno_name = "EasilyBurn", assoc="", maf = 0.01, hwe = 1e-4, geno=0.02)
hubcollr = collrinfo()
bedcollinfo(hubcollr, genbed=FALSE)
taskinfo(hubcollr, taskName, plinkArgs, initGwasArgs, TRUE, "bin")
taskBedsPlinkOut(hubcollr, taskName, hubcollr[[taskName]]$fullGwasOut, 1e-3, 1, 200)
taskAnalyze(hubcollr, taskName)
readcoll.task(hubcollr, taskName)
globalMin = min(hubcollr[[taskName]]$minPvals)
if(globalMin > 1) globalMin = 1
if(globalMin < 1e-300) globalMin = 1e-300
globalMinVec = c(globalMinVec, globalMin)
message("*****************************************************")
print(globalMinVec)
message("*****************************************************")
unlink(hubcollr[[taskName]]$taskPath, recursive = TRUE, force = TRUE)
}
perm
# erf4 huid07, winszie=200
setwd("/media/data1/kaiyin/erf4/")
require(collr2)
plinkArgs = getPlinkParam(allow_no_sex = "", pheno = "erf4.pheno/erf_all_phdata.txt", pheno_name = "huid07", covar = "erf4.pheno/erf_all_phdata.txt", covar_name = "sex,age", linear = "hide-covar")
initGwasArgs = getPlinkParam(maf=.01, hwe=1e-4, allow_no_sex = "", pheno = "erf4.pheno/erf_all_phdata.txt", pheno_name = "huid07", assoc= "")
crt1(".", "huid07b", plinkArgs, initGwasArgs, TRUE, "con", 0.001, 200)
shuffleLines = function(phenoFileOrig, phenoFile, shuffleVar) {
tab = read.table(phenoFileOrig, head=TRUE)
# print(head(tab, 10))
tab[, shuffleVar] = sample(tab[, shuffleVar])
# print(head(tab, 10))
write.table(tab, phenoFile, quote = FALSE, row.names = FALSE)
}
permute = function(plinkArgs, initGwasArgs, phenoFileOrig, phenoFileShuffle = "~/Desktop/shuffle_pheno.csv", wDir=".", n=100) {
setwd(wDir)
globalMinVec = numeric(0)
plinkArgs$pheno = phenoFileShuffle
if(!is.null(plinkArgs$covar)) plinkArgs$covar = phenoFileShuffle
initGwasArgs$pheno = phenoFileShuffle
taskName = "permutation_analysis"
for(i in 1:n) {
shuffleLines(phenoFileOrig, phenoFileShuffle, plinkArgs$pheno_name)
hubcollr = collrinfo()
bedcollinfo(hubcollr, genbed=FALSE)
taskinfo(hubcollr, taskName, plinkArgs, initGwasArgs, TRUE)
taskBedsPlinkOut(hubcollr, taskName, hubcollr[[taskName]]$fullGwasOut, 1e-3, 1, 200)
taskAnalyze(hubcollr, taskName)
readcoll.task(hubcollr, taskName)
globalMin = min(hubcollr[[taskName]]$minPvals)
if(globalMin > 1) globalMin = 1
if(globalMin < 1e-300) globalMin = 1e-300
globalMinVec = c(globalMinVec, globalMin)
message("*****************************************************")
print(globalMinVec)
message("*****************************************************")
unlink(hubcollr[[taskName]]$taskPath, recursive = TRUE, force = TRUE)
}
globalMinVec
}
require(collr2)
setwd("/media/data1/kaiyin/erf4")
phenoFileOrig = "erf4.pheno/erf_all_phdata.txt"
require(collr2)
plinkArgs = getPlinkParam(allow_no_sex = "", pheno = "erf4.pheno/erf_all_phdata.txt", pheno_name = "huid07", covar = "erf4.pheno/erf_all_phdata.txt", covar_name = "sex,age", linear = "hide-covar")
initGwasArgs = getPlinkParam(maf=.01, hwe=1e-4, allow_no_sex = "", pheno = "erf4.pheno/erf_all_phdata.txt", pheno_name = "huid07", assoc= "")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.