-
Notifications
You must be signed in to change notification settings - Fork 0
/
run_t1d.r
72 lines (55 loc) · 2.28 KB
/
run_t1d.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
options(mc.cores=4)
options(boot.ncpus=4)
options(boot.parallel="snow")
#cl<-makeCluster(getOption("mc.cores", 2L)) #number of CPU cores
#registerDoSNOW(cl)
#power.nistal()
#proc.choc()
#meta = load.meta.t1d("annotation20130819.csv")
#mgrast.dir = "../BATCH_01_02_META/BATCH_01-02_METAGENOMICS_MGRAST"
#mgr = read.mgrast.summary(paste(mgrast.dir,"cog.tsv",sep="/"),file_name.id.map=paste(mgrast.dir,"mgrast_to_samp_id.tsv",sep="/"))
#mgr.cnt = mgrast.to.abund.df(mgr,"level.2")
#taxa.meta = read.t1d(3)
#taxa.meta.data = taxa.meta$data
#taxa.meta.attr.names = taxa.meta$attr.names
#sink("analysis.log",split=T)
# panderOptions("table.style","rmarkdown")
# panderOptions("table.split.table",180)
# panderOptions("table.alignment.default","left")
# panderOptions("evals.messages",F)
#
# evalsOptions("cache",F)
# evalsOptions("cache.mode","environment")
# #evalsOptions("output",c("result"))
# evalsOptions("output",c("all"))
# evalsOptions("graph.unify",F)
# evalsOptions("res",100)
# evalsOptions("hi.res",T)
MGSAT_SRC = "~/work/mgsat"
source(paste(MGSAT_SRC,"dependencies.r",sep="/"),local=T)
#install_required_packages()
load_required_packages()
source(paste(MGSAT_SRC,"report_pandoc.r",sep="/"),local=T)
source(paste(MGSAT_SRC,"power_and_tests.r",sep="/"),local=T)
set_trace_options(try.debug=F)
setwd("~/work/T1D/2014-10-09")
report <- PandocAT$new("atovtchi@jcvi.org","T1D 16S Analysis MiSeq V1-V3")
proc.t1d()
#proc.t1d.mg()
#print(report$body)
report$save("report.md")
Pandoc.convert("report.md",format="html",open=F)
Pandoc.convert("report.md",format="docx",open=F)
#Pandoc.convert("report.md",format="pdf",open=F)
#print(plot.stability.selection.c060.at(s),rank="mean")
#spca.res = spca(count,K=6,type="predictor",sparse="varnum",trace=T,para=c(7,4,4,1,1,1))
##This has selected a single variable (Gordonibacter) from 16S level 3 with normalization "ident".
##Note that by default this function does not scale the predictor variables.
#hc.res <- get.biom(X = m_a$count, Y = m_a$attr$T1D, fmethod = c("studentt", "pls", "vip"), type = "HC")
#taxa.meta$data[taxa.meta$data$Batch==3,c("Gordonibacter_0.1.1.1.3.1.6","T1D")]
#taxa.meta$data[taxa.meta$data$Batch==3,c("Akkermansia_0.1.8.2.1.1.1","T1D")]
#taxa.meta = read.mr_oralc(3)
#proc.mr_oralc()
#power.pediatric.cancer.2013()
#stopCluster(cl)
#sink(NULL)