Mercurial > repos > george-weingart > maaslin

comparison maaslin-4450aa4ecc84/src/test-BoostGLM/test-BoostGLM.R @ 1:a87d5a5f2776

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded the version running on the prod server
author george-weingart
date Sun, 08 Feb 2015 23:08:38 -0500
parents
children
comparison
equal deleted inserted replaced
0:e0b5980139d9 1:a87d5a5f2776
1 c_strDir <- file.path(getwd( ),"..")
2
3 source(file.path(c_strDir,"lib","Constants.R"))
4 source(file.path(c_strDir,"lib","Utility.R"))
5 source(file.path(c_strDir,"lib","AnalysisModules.R"))
6
7 # General setup
8 covX1 = c(44.4, 45.9, 41.9, 53.3, 44.7, 44.1, 50.7, 45.2, 60.1)
9 covX2 = c(144.4, 245.9, 141.9, 253.3, 144.7, 244.1, 150.7, 245.2, 160.1)
10 covX3 = as.factor(c(1,2,3,1,2,3,1,2,3))
11 covX4 = as.factor(c(1,1,1,1,2,2,2,2,2))
12 covX5 = as.factor(c(1,2,1,2,1,2,1,2,1))
13 covY = c(.26, .31, .25, .50, .36, .40, .52, .28, .38)
14 frmeTmp = data.frame(Covariate1=covX1, Covariate2=covX2, Covariate3=covX3, Covariate4=covX4, Covariate5=covX5, adCur= covY)
15 iTaxon = 6
16
17 lsCov1 = list()
18 lsCov1$name = "Covariate1"
19 lsCov1$orig = "Covariate1"
20 lsCov1$taxon = "adCur"
21 lsCov1$data = covY
22 lsCov1$factors = "Covariate1"
23 lsCov1$metadata = covX1
24 vdCoef = c()
25 vdCoef["(Intercept)"]=round(0.0345077486,5)
26 vdCoef["Covariate1"]= round(0.0052097355,5)
27 vdCoef["Covariate2"]= round(0.0005806568,5)
28 vdCoef["Covariate32"]=round(-0.1333421874,5)
29 vdCoef["Covariate33"]=round(-0.1072006419,5)
30 vdCoef["Covariate42"]=round(0.0849198280,5)
31 lsCov1$value = c(Covariate1=round(0.005209736,5))
32 lsCov1$std = round(0.0063781728,5)
33 lsCov1$allCoefs = vdCoef
34 lsCov2 = list()
35 lsCov2$name = "Covariate2"
36 lsCov2$orig = "Covariate2"
37 lsCov2$taxon = "adCur"
38 lsCov2$data = covY
39 lsCov2$factors = "Covariate2"
40 lsCov2$metadata = covX2
41 lsCov2$value = c(Covariate2=round(0.0005806568,5))
42 lsCov2$std = round(0.0006598436,5)
43 lsCov2$allCoefs = vdCoef
44 lsCov3 = list()
45 lsCov3$name = "Covariate3"
46 lsCov3$orig = "Covariate32"
47 lsCov3$taxon = "adCur"
48 lsCov3$data = covY
49 lsCov3$factors = "Covariate3"
50 lsCov3$metadata = covX3
51 lsCov3$value = c(Covariate32=round(-0.1333422,5))
52 lsCov3$std = round(0.0895657826,5)
53 lsCov3$allCoefs = vdCoef
54 lsCov4 = list()
55 lsCov4$name = "Covariate3"
56 lsCov4$orig = "Covariate33"
57 lsCov4$taxon = "adCur"
58 lsCov4$data = covY
59 lsCov4$factors = "Covariate3"
60 lsCov4$metadata = covX3
61 lsCov4$value = c(Covariate33=round(-0.1072006,5))
62 lsCov4$std = round(0.0792209541,5)
63 lsCov4$allCoefs = vdCoef
64 lsCov5 = list()
65 lsCov5$name = "Covariate4"
66 lsCov5$orig = "Covariate42"
67 lsCov5$taxon = "adCur"
68 lsCov5$data = covY
69 lsCov5$factors = "Covariate4"
70 lsCov5$metadata = covX4
71 lsCov5$value = c(Covariate42=round(0.08491983,5))
72 lsCov5$std = round(0.0701018621,5)
73 lsCov5$allCoefs = vdCoef
74
75 context("Test funcClean")
76
77 context("Test funcBugHybrid")
78 # multiple covariates, one call lm
79 aiMetadata = c(1:5)
80 aiData = c(iTaxon)
81 dFreq = 0.5 / length( aiMetadata )
82 dSig = 0.25
83 dMinSamp = 0.1
84 adP = c()
85 lsSig = list()
86 funcReg = NA
87 funcAnalysis = funcLM
88 funcGetResult = funcGetLMResults
89 lsData = list(frmeData=frmeTmp, aiMetadata=aiMetadata, aiData=aiData, lsQCCounts=list())
90 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
91
92 adPExpected = round(c(0.4738687,0.4436566,0.4665972,0.5378693,0.3124672),5)
93 QCExpected = list(iLms=numeric(0))
94 lsSigExpected = list()
95 lsSigExpected[[1]] = lsCov1
96 lsSigExpected[[2]] = lsCov2
97 lsSigExpected[[3]] = lsCov3
98 lsSigExpected[[4]] = lsCov4
99 lsSigExpected[[5]] = lsCov5
100 expectedReturn = list(adP=adPExpected,lsSig=lsSigExpected,lsQCCounts=QCExpected)
101 receivedReturn = funcBugHybrid(iTaxon=iTaxon,frmeData=frmeTmp,lsData=lsData,aiMetadata=aiMetadata,dFreq=dFreq,dSig=dSig,dMinSamp=dMinSamp,adP=adP,lsSig=lsSig, strLog=NA,funcReg=funcReg,lsNonPenalizedPredictors=NULL,funcAnalysis=funcAnalysis,lsRandomCovariates=NULL,funcGetResult=funcGetResult)
102 receivedReturn$adP = round(receivedReturn$adP,5)
103
104 vCoefs=receivedReturn$lsSig[[1]]$allCoefs
105 vCoefs[1]=round(vCoefs[1],5)
106 vCoefs[2]=round(vCoefs[2],5)
107 vCoefs[3]=round(vCoefs[3],5)
108 vCoefs[4]=round(vCoefs[4],5)
109 vCoefs[5]=round(vCoefs[5],5)
110 vCoefs[6]=round(vCoefs[6],5)
111 receivedReturn$lsSig[[1]]$allCoefs=vCoefs
112 receivedReturn$lsSig[[2]]$allCoefs=vCoefs
113 receivedReturn$lsSig[[3]]$allCoefs=vCoefs
114 receivedReturn$lsSig[[4]]$allCoefs=vCoefs
115 receivedReturn$lsSig[[5]]$allCoefs=vCoefs
116 vValue=c()
117 vValue[receivedReturn$lsSig[[1]]$orig]=round(receivedReturn$lsSig[[1]]$value[[1]],5)
118 receivedReturn$lsSig[[1]]$value=vValue
119 vValue=c()
120 vValue[receivedReturn$lsSig[[2]]$orig]=round(receivedReturn$lsSig[[2]]$value[[1]],5)
121 receivedReturn$lsSig[[2]]$value=vValue
122 vValue=c()
123 vValue[receivedReturn$lsSig[[3]]$orig]=round(receivedReturn$lsSig[[3]]$value[[1]],5)
124 receivedReturn$lsSig[[3]]$value=vValue
125 vValue=c()
126 vValue[receivedReturn$lsSig[[4]]$orig]=round(receivedReturn$lsSig[[4]]$value[[1]],5)
127 receivedReturn$lsSig[[4]]$value=vValue
128 vValue=c()
129 vValue[receivedReturn$lsSig[[5]]$orig]=round(receivedReturn$lsSig[[5]]$value[[1]],5)
130 receivedReturn$lsSig[[5]]$value=vValue
131 receivedReturn$lsSig[[1]]$std=round(receivedReturn$lsSig[[1]]$std,5)
132 receivedReturn$lsSig[[2]]$std=round(receivedReturn$lsSig[[2]]$std,5)
133 receivedReturn$lsSig[[3]]$std=round(receivedReturn$lsSig[[3]]$std,5)
134 receivedReturn$lsSig[[4]]$std=round(receivedReturn$lsSig[[4]]$std,5)
135 receivedReturn$lsSig[[5]]$std=round(receivedReturn$lsSig[[5]]$std,5)
136 test_that("funcBugHybrid works with the lm option with multiple covariates.",{expect_equal(receivedReturn,expectedReturn)})
137
138
139 # single covariate, single call lm
140 aiMetadata = c(1)
141 dFreq = 0.5 / length( aiMetadata )
142 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
143 adPExpected = round(c(0.1081731),5)
144 QCExpected = list(iLms=numeric(0))
145 lsSigExpected = list()
146 lsSigExpected[[1]] = lsCov1
147 lsSigExpected[[1]]$std=round(0.005278468,5)
148 vdCoef = c()
149 vdCoef["(Intercept)"]=round(-0.102410716,5)
150 vdCoef["Covariate1"]= round(0.009718095,5)
151 lsSigExpected[[1]]$allCoefs= vdCoef
152 lsSigExpected[[1]]$value = c(Covariate1=round(0.009718095,5))
153
154 expectedReturn = list(adP=adPExpected,lsSig=lsSigExpected,lsQCCounts=QCExpected)
155 receivedReturn = funcBugHybrid(iTaxon=iTaxon,frmeData=frmeTmp,lsData=lsData,aiMetadata=aiMetadata,dFreq=dFreq,dSig=dSig,dMinSamp=dMinSamp,adP=adP,lsSig=lsSig, strLog=NA,funcReg=funcReg,lsNonPenalizedPredictors=NULL,funcAnalysis=funcAnalysis,lsRandomCovariates=NULL,funcGetResult=funcGetResult)
156 receivedReturn$adP = round(receivedReturn$adP,5)
157
158 vCoefs=receivedReturn$lsSig[[1]]$allCoefs
159 vCoefs[1]=round(vCoefs[1],5)
160 vCoefs[2]=round(vCoefs[2],5)
161 receivedReturn$lsSig[[1]]$allCoefs=vCoefs
162 vValue=c()
163 vValue[receivedReturn$lsSig[[1]]$orig]=round(receivedReturn$lsSig[[1]]$value[[1]],5)
164 receivedReturn$lsSig[[1]]$value=vValue
165 receivedReturn$lsSig[[1]]$std=round(0.005278468,5)
166 test_that("funcBugHybrid works with the lm option with 1 covariates.",{expect_equal(receivedReturn,expectedReturn)})
167
168
169 # multiple covariate, single call univariate
170 funcReg = NA
171 funcAnalysis = funcDoUnivariate
172 funcGetResult = NA
173 aiMetadata = c(3,1,2)
174 dFreq = 0.5 / length( aiMetadata )
175 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
176 adPExpected = round(c(1.0,1.0,0.09679784,0.21252205),5)
177 QCExpected = list(iLms=numeric(0))
178 lsSigExpected = list()
179 lsCov1 = list()
180 lsCov1$name = "Covariate3"
181 lsCov1$orig = "Covariate32"
182 lsCov1$taxon = "adCur"
183 lsCov1$data = covY
184 lsCov1$factors = "Covariate3"
185 lsCov1$metadata = frmeTmp[["Covariate3"]]
186 vdCoef = c(Covariate32=NA)
187 lsCov1$value = vdCoef
188 lsCov1$std = sd(frmeTmp[["Covariate3"]])
189 lsCov1$allCoefs = vdCoef
190 lsCov2 = list()
191 lsCov2$name = "Covariate3"
192 lsCov2$orig = "Covariate33"
193 lsCov2$taxon = "adCur"
194 lsCov2$data = covY
195 lsCov2$factors = "Covariate3"
196 lsCov2$metadata = frmeTmp[["Covariate3"]]
197 vdCoef = c(Covariate33=NA)
198 lsCov2$value = vdCoef
199 lsCov2$std = sd(frmeTmp[["Covariate3"]])
200 lsCov2$allCoefs = vdCoef
201 lsCov3 = list()
202 lsCov3$name = "Covariate1"
203 lsCov3$orig = "Covariate1"
204 lsCov3$taxon = "adCur"
205 lsCov3$data = covY
206 lsCov3$factors = "Covariate1"
207 lsCov3$metadata = frmeTmp[["Covariate1"]]
208 vdCoef = c(Covariate1=0.6)
209 lsCov3$value = vdCoef
210 lsCov3$std = sd(frmeTmp[["Covariate1"]])
211 lsCov3$allCoefs = vdCoef
212 lsCov4 = list()
213 lsCov4$name = "Covariate2"
214 lsCov4$orig = "Covariate2"
215 lsCov4$taxon = "adCur"
216 lsCov4$data = covY
217 lsCov4$factors = "Covariate2"
218 lsCov4$metadata = frmeTmp[["Covariate2"]]
219 vdCoef = c(Covariate2=0.46666667)
220 lsCov4$value = vdCoef
221 lsCov4$std = sd(frmeTmp[["Covariate2"]])
222 lsCov4$allCoefs = vdCoef
223
224 lsSigExpected = list()
225 lsSigExpected[[1]] = lsCov1
226 lsSigExpected[[2]] = lsCov2
227 lsSigExpected[[3]] = lsCov3
228 lsSigExpected[[4]] = lsCov4
229
230 expectedReturn = list(adP=adPExpected,lsSig=lsSigExpected,lsQCCounts=QCExpected)
231 receivedReturn = funcBugHybrid(iTaxon=iTaxon,frmeData=frmeTmp,lsData=lsData,aiMetadata=aiMetadata,dFreq=dFreq,dSig=dSig,dMinSamp=dMinSamp,adP=adP,lsSig=lsSig, strLog=NA,funcReg=funcReg,lsNonPenalizedPredictors=NULL,funcAnalysis=funcAnalysis,lsRandomCovariates=NULL,funcGetResult=funcGetResult)
232 receivedReturn$adP = round(receivedReturn$adP,5)
233 test_that("funcBugHybrid works with the univariate option with 3 covariates.",{expect_equal(receivedReturn,expectedReturn)})
234
235
236 # single covariate, single call univariate
237 funcReg = NA
238 funcAnalysis = funcDoUnivariate
239 funcGetResult = NA
240 aiMetadata = c(1)
241 dFreq = 0.5 / length( aiMetadata )
242 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
243 adPExpected = round(c(0.09679784),5)
244 QCExpected = list(iLms=numeric(0))
245 lsSigExpected = list()
246 lsSigExpected[[1]] = lsCov3
247
248 expectedReturn = list(adP=adPExpected,lsSig=lsSigExpected,lsQCCounts=QCExpected)
249 receivedReturn = funcBugHybrid(iTaxon=iTaxon,frmeData=frmeTmp,lsData=lsData,aiMetadata=aiMetadata,dFreq=dFreq,dSig=dSig,dMinSamp=dMinSamp,adP=adP,lsSig=lsSig, strLog=NA,funcReg=funcReg,lsNonPenalizedPredictors=NULL,funcAnalysis=funcAnalysis,lsRandomCovariates=NULL,funcGetResult=funcGetResult)
250 receivedReturn$adP = round(receivedReturn$adP,5)
251 test_that("funcBugHybrid works with the univariate option with 1 covariates.",{expect_equal(receivedReturn,expectedReturn)})
252
253
254 context("Test funcBugs")
255 #One LM run
256 frmeData=frmeTmp
257 aiMetadata=c(1)
258 aiData=c(iTaxon)
259 strData=NA
260 dFreq= 0.5 / length( aiMetadata )
261 dSig=0.25
262 dMinSamp=0.1
263 strDirOut=NA
264 funcReg=NA
265 lsNonPenalizedPredictors=NULL
266 lsRandomCovariates=NULL
267 funcAnalysis=funcLM
268 funcGetResults=funcGetLMResults
269 fDoRPlot=FALSE
270 lsData = list(frmeData=frmeData, aiMetadata=aiMetadata, aiData=aiData, lsQCCounts=list())
271 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
272 QCExpected = list(iLms=numeric(0))
273
274 expectedReturn = list(aiReturnBugs=aiData,lsQCCounts=QCExpected)
275 receivedReturn = funcBugs(frmeData=frmeData, lsData=lsData, aiMetadata=aiMetadata, aiData=aiData, strData=strData, dFreq=dFreq, dSig=dSig, dMinSamp=dMinSamp,strDirOut=strDirOut, funcReg=funcReg,lsNonPenalizedPredictors=lsNonPenalizedPredictors,funcAnalysis=funcAnalysis,lsRandomCovariates=lsRandomCovariates,funcGetResults=funcGetResults,fDoRPlot=fDoRPlot)
276
277 test_that("funcBugs works with the lm option with 1 covariate.",{expect_equal(receivedReturn,expectedReturn)})
278
279 #multiple LM run
280 frmeData=frmeTmp
281 aiMetadata=c(1:5)
282 aiData=c(iTaxon)
283 strData=NA
284 dFreq= 0.5 / length( aiMetadata )
285 dSig=0.25
286 dMinSamp=0.1
287 strDirOut=NA
288 funcReg=NA
289 lsNonPenalizedPredictors=NULL
290 lsRandomCovariates=NULL
291 funcAnalysis=funcLM
292 funcGetResults=funcGetLMResults
293 fDoRPlot=FALSE
294 lsData = list(frmeData=frmeData, aiMetadata=aiMetadata, aiData=aiData, lsQCCounts=list())
295 lsData$astrMetadata = names(frmeTmp)[aiMetadata]
296 QCExpected = list(iLms=numeric(0))
297
298 expectedReturn = list(aiReturnBugs=aiData,lsQCCounts=QCExpected)
299 receivedReturn = funcBugs(frmeData=frmeData, lsData=lsData, aiMetadata=aiMetadata, aiData=aiData, strData=strData, dFreq=dFreq, dSig=dSig, dMinSamp=dMinSamp,strDirOut=strDirOut, funcReg=funcReg,lsNonPenalizedPredictors=lsNonPenalizedPredictors,funcAnalysis=funcAnalysis,lsRandomCovariates=lsRandomCovariates,funcGetResults=funcGetResults,fDoRPlot=fDoRPlot)
300
301 print("START START")
302 print(expectedReturn)
303 print("RECEIVED")
304 print(receivedReturn)
305 print("STOP STOP")
306
307 test_that("funcBugs works with the lm option with multiple covariates.",{expect_equal(receivedReturn,expectedReturn)})