processedTongue = processDataCube(parafac4microbiome::vanderPloeg2024$tongue, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedLowling = processDataCube(parafac4microbiome::vanderPloeg2024$lower_jaw_lingual, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedLowinter = processDataCube(parafac4microbiome::vanderPloeg2024$lower_jaw_interproximal, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedUpling = processDataCube(parafac4microbiome::vanderPloeg2024$upper_jaw_lingual, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedUpinter = processDataCube(parafac4microbiome::vanderPloeg2024$upper_jaw_interproximal, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedSaliva = processDataCube(parafac4microbiome::vanderPloeg2024$saliva, sparsityThreshold=0.50, considerGroups=TRUE, groupVariable="RFgroup", CLR=TRUE, centerMode=1, scaleMode=2)
processedMetabolomics = parafac4microbiome::vanderPloeg2024$metabolomics
assessment_tongue = assessModelQuality(processedTongue$data, numRepetitions=10, numCores=10)
assessment_tongue$plots$overview
assessment_lowling = parafac4microbiome::assessModelQuality(processedLowling$data, numRepetitions=10, numCores=10)
assessment_lowling$plots$overview
assessment_lowinter = parafac4microbiome::assessModelQuality(processedLowinter$data, numRepetitions=10, numCores=10)
assessment_lowinter$plots$overview
assessment_upling = parafac4microbiome::assessModelQuality(processedUpling$data, numRepetitions=10, numCores=10)
assessment_upling$plots$overview
assessment_upinter = parafac4microbiome::assessModelQuality(processedUpinter$data, numRepetitions=10, numCores=10)
assessment_upinter$plots$overview
assessment_saliva = parafac4microbiome::assessModelQuality(processedSaliva$data, numRepetitions=10, numCores=10)
assessment_saliva$plots$overview
assessment_metab = parafac4microbiome::assessModelQuality(processedMetabolomics$data, numRepetitions=10, numCores=10)
assessment_metab$plots$overview
colourCols = c("RFgroup", "Phylum", "")
legendTitles = c("RF group", "Phylum", "")
xLabels = c("Subject index", "Feature index", "Time index")
legendColNums = c(3,5,0)
arrangeModes = c(TRUE, TRUE, FALSE)
continuousModes = c(FALSE,FALSE,TRUE)
stability_tongue = assessModelStability(processedTongue, maxNumComponents=3, numFolds=10, colourCols=colourCols, legendTitles=legendTitles, xLabels = xLabels, legendColNums=legendColNums, arrangeModes=arrangeModes, numCores=parallel::detectCores())
stability_tongue$modelPlots[[2]]
stability_tongue$modelPlots[[3]]
stability_upinter = assessModelStability(processedUpinter, maxNumComponents=3, numFolds=10, colourCols=colourCols, legendTitles=legendTitles, xLabels = xLabels, legendColNums=legendColNums, arrangeModes=arrangeModes, numCores=parallel::detectCores())
stability_upinter$modelPlots[[2]]
stability_upinter$modelPlots[[3]]
cp_tongue = parafac4microbiome::parafac(processedTongue$data, nfac=2, nstart=100)
cp_lowling = parafac4microbiome::parafac(processedLowling$data, nfac=1, nstart=100)
cp_lowinter = parafac4microbiome::parafac(processedLowinter$data, nfac=2, nstart=100)
cp_upling = parafac4microbiome::parafac(processedUpling$data, nfac=2, nstart=100)
cp_upinter = parafac4microbiome::parafac(processedUpinter$data, nfac=2, nstart=100)
cp_saliva = parafac4microbiome::parafac(processedSaliva$data, nfac=2, nstart=100)
cp_metabolomics = parafac4microbiome::parafac(processedMetabolomics$data, nfac=1, nstart=100)
testMetadata = function(model, comp, metadata){
transformedSubjectLoadings = transformPARAFACloadings(model$Fac, 1)[,comp]
transformedSubjectLoadings = transformedSubjectLoadings / norm(transformedSubjectLoadings, "2")
metadata = metadata$mode1 %>% left_join(vanderPloeg2024$red_fluorescence %>% select(subject,day,Area_delta_R30,plaquepercent,bomppercent) %>% filter(day==14), by="subject")
result = lm(transformedSubjectLoadings ~ plaquepercent + bomppercent + Area_delta_R30 + gender + age, data=metadata)
# Extract coefficients and confidence intervals
coef_estimates <- summary(result)$coefficients
conf_intervals <- confint(result)
# Remove intercept
coef_estimates <- coef_estimates[rownames(coef_estimates) != "(Intercept)", ]
conf_intervals <- conf_intervals[rownames(conf_intervals) != "(Intercept)", ]
# Combine into a clean data frame
summary_table <- data.frame(
Term = rownames(coef_estimates),
Estimate = coef_estimates[, "Estimate"] * 1e3,
CI = paste0(
conf_intervals[, 1], " – ",
conf_intervals[, 2]
),
P_value = coef_estimates[, "Pr(>|t|)"],
P_adjust = p.adjust(coef_estimates[, "Pr(>|t|)"], "BH"),
row.names = NULL
)
return(summary_table)
}
testMetadata(cp_tongue, 1, processedTongue)
#> Term Estimate CI
#> 1 plaquepercent -7.5822683 -0.0128447981932743 – -0.00231973832740774
#> 2 bomppercent 0.9231415 -0.00236470070402032 – 0.00421098379021022
#> 3 Area_delta_R30 1.8514918 -0.00477870837570338 – 0.0084816919298291
#> 4 gender 4.8154398 -0.0885416859076374 – 0.0981725654207313
#> 5 age 14.8619327 0.0065466568449098 – 0.0231772086332818
#> P_value P_adjust
#> 1 0.0060091824 0.015022956
#> 2 0.5723172847 0.717991555
#> 3 0.5743932438 0.717991555
#> 4 0.9172000037 0.917200004
#> 5 0.0009008391 0.004504195
testMetadata(cp_tongue, 2, processedTongue)
#> Term Estimate CI
#> 1 plaquepercent -1.470831 -0.00774087652102776 – 0.00479921503852629
#> 2 bomppercent -2.022641 -0.00593994325341581 – 0.00189466164910002
#> 3 Area_delta_R30 -4.994439 -0.0128939967625232 – 0.00290511849828337
#> 4 gender -42.654012 -0.153884438874493 – 0.0685764152176233
#> 5 age 3.393808 -0.0065134346803709 – 0.0133010504408422
#> P_value P_adjust
#> 1 0.6368744 0.6368744
#> 2 0.3017239 0.6142220
#> 3 0.2077427 0.6142220
#> 4 0.4415093 0.6142220
#> 5 0.4913776 0.6142220
testMetadata(cp_lowling, 1, processedLowling)
#> Term Estimate CI
#> 1 plaquepercent -4.045863 -0.00898804283492606 – 0.00089631674972345
#> 2 bomppercent -2.596336 -0.00568403462309891 – 0.00049136303586988
#> 3 Area_delta_R30 -10.214824 -0.0164414185575453 – -0.00398822907200625
#> 4 gender -82.762670 -0.170436793375464 – 0.00491145387711936
#> 5 age -2.443801 -0.0102528940991745 – 0.00536529297646459
#> P_value P_adjust
#> 1 0.105459929 0.13182491
#> 2 0.096672023 0.13182491
#> 3 0.002053051 0.01026526
#> 4 0.063511447 0.13182491
#> 5 0.529359622 0.52935962
testMetadata(cp_lowinter, 1, processedLowinter)
#> Term Estimate CI
#> 1 plaquepercent -2.434126 -0.00846635949340225 – 0.00359810745608451
#> 2 bomppercent 3.424065 -0.000344660250091408 – 0.00719279124775479
#> 3 Area_delta_R30 -7.306946 -0.0149068863574801 – 0.000292994998243587
#> 4 gender -78.160031 -0.185171674697412 – 0.0288516128522022
#> 5 age -5.976242 -0.015507719707431 – 0.00355523527425291
#> P_value P_adjust
#> 1 0.41822203 0.4182220
#> 2 0.07359514 0.1839878
#> 3 0.05899468 0.1839878
#> 4 0.14708562 0.2451427
#> 5 0.21145149 0.2643144
testMetadata(cp_lowinter, 2, processedLowinter)
#> Term Estimate CI
#> 1 plaquepercent 2.8311345 -0.00314110864245655 – 0.00880337766140032
#> 2 bomppercent 3.9200707 0.000188824738963201 – 0.00765131658041496
#> 3 Area_delta_R30 0.3307337 -0.00719362590923236 – 0.0078550932288655
#> 4 gender 19.8248632 -0.0861225539988762 – 0.125772280439005
#> 5 age 3.2431239 -0.0061935634485384 – 0.0126798112014249
#> P_value P_adjust
#> 1 0.34246585 0.8166329
#> 2 0.04002742 0.2001371
#> 3 0.92940501 0.9294050
#> 4 0.70633476 0.8829185
#> 5 0.48997974 0.8166329
testMetadata(cp_upling, 1, processedUpling)
#> Term Estimate CI P_value
#> 1 plaquepercent 2.511420 -0.00356746805080028 – 0.00859030851362675 0.4073212
#> 2 bomppercent 3.020891 -0.00077698350686279 – 0.00681876453237555 0.1153387
#> 3 Area_delta_R30 3.798073 -0.0038606477132654 – 0.0114567933392421 0.3209597
#> 4 gender 34.561948 -0.0732773509874726 – 0.14240124637675 0.5195281
#> 5 age -0.886911 -0.0104921073081138 – 0.00871828535120744 0.8523878
#> P_adjust
#> 1 0.6494101
#> 2 0.5766934
#> 3 0.6494101
#> 4 0.6494101
#> 5 0.8523878
testMetadata(cp_upling, 2, processedUpling)
#> Term Estimate CI
#> 1 plaquepercent 2.5306859 -0.00356824415515762 – 0.0086296158954944
#> 2 bomppercent 0.2034754 -0.0036069199493453 – 0.00401387083190149
#> 3 Area_delta_R30 8.1783956 0.00049442470478565 – 0.0158623664746058
#> 4 gender 84.5889919 -0.0236058467630009 – 0.192783830463618
#> 5 age 0.6422043 -0.00899465978137149 – 0.010279068431573
#> P_value P_adjust
#> 1 0.40529959 0.6754993
#> 2 0.91429122 0.9142912
#> 3 0.03764326 0.1882163
#> 4 0.12146608 0.3036652
#> 5 0.89315990 0.9142912
testMetadata(cp_upinter, 1, processedUpinter)
#> Term Estimate CI
#> 1 plaquepercent -2.814274 -0.00826061185340105 – 0.00263206420300627
#> 2 bomppercent -3.192434 -0.00659511342869162 – 0.000210244636499445
#> 3 Area_delta_R30 -7.022575 -0.0138843530784991 – -0.000160797417321151
#> 4 gender -110.454309 -0.207072184351565 – -0.0138364326547657
#> 5 age -1.389587 -0.00999529601516111 – 0.0072161227870755
#> P_value P_adjust
#> 1 0.30136208 0.3767026
#> 2 0.06506720 0.1084453
#> 3 0.04513921 0.1084453
#> 4 0.02624554 0.1084453
#> 5 0.74501059 0.7450106
testMetadata(cp_upinter, 2, processedUpinter)
#> Term Estimate CI P_value
#> 1 plaquepercent -2.710847 -0.00896645118894094 – 0.00354475799069917 0.3850003
#> 2 bomppercent 2.584954 -0.00132332573260628 – 0.00649323448226151 0.1880057
#> 3 Area_delta_R30 -1.600749 -0.00948211239065241 – 0.00628061429618728 0.6826150
#> 4 gender 21.046166 -0.0899280749044152 – 0.132020406620614 0.7025612
#> 5 age -7.790367 -0.0176747911938898 – 0.00209405714035475 0.1185832
#> P_adjust
#> 1 0.6416672
#> 2 0.4700141
#> 3 0.7025612
#> 4 0.7025612
#> 5 0.4700141
testMetadata(cp_saliva, 1, processedSaliva)
#> Term Estimate CI
#> 1 plaquepercent -3.897532 -0.00981027366469069 – 0.00201521022558586
#> 2 bomppercent 2.211191 -0.0014828807091094 – 0.00590526264067308
#> 3 Area_delta_R30 3.876556 -0.00357283918714432 – 0.0113259501964904
#> 4 gender 69.061034 -0.0358308334679268 – 0.173952901665905
#> 5 age 10.372091 0.00102942100161564 – 0.0197147610154618
#> P_value P_adjust
#> 1 0.18946185 0.2905401
#> 2 0.23243207 0.2905401
#> 3 0.29800575 0.2980057
#> 4 0.18996862 0.2905401
#> 5 0.03057498 0.1528749
testMetadata(cp_saliva, 2, processedSaliva)
#> Term Estimate CI
#> 1 plaquepercent -4.132439 -0.00974061125477911 – 0.00147573233866911
#> 2 bomppercent -1.461212 -0.00496499904532008 – 0.00204257502925521
#> 3 Area_delta_R30 -6.286724 -0.0133523944344635 – 0.000778946185827744
#> 4 gender -73.637107 -0.173125908829778 – 0.0258516954196962
#> 5 age 10.065789 0.00120436754605748 – 0.0189272104927239
#> P_value P_adjust
#> 1 0.14364164 0.1795520
#> 2 0.40295215 0.4029522
#> 3 0.07947659 0.1795520
#> 4 0.14191360 0.1795520
#> 5 0.02715158 0.1357579
testMetadata(cp_metabolomics, 1, processedMetabolomics)
#> Term Estimate CI P_value
#> 1 plaquepercent -4.322341 -0.0104408654355027 – 0.00179618427456725 0.1602392
#> 2 bomppercent -1.786318 -0.00552918075883336 – 0.00195654395477168 0.3389415
#> 3 Area_delta_R30 -3.987202 -0.0115380347230205 – 0.00356363052776998 0.2907738
#> 4 gender 34.644787 -0.0718282432667193 – 0.141117817257719 0.5129001
#> 5 age 7.326636 -0.00213977971574463 – 0.0167930523653612 0.1250078
#> P_adjust
#> 1 0.4005979
#> 2 0.4236768
#> 3 0.4236768
#> 4 0.5129001
#> 5 0.4005979
colourCols = c("RFgroup", "Phylum", "")
legendTitles = c("RF group", "Phylum", "")
xLabels = c("Subject index", "Feature index", "Time index")
legendColNums = c(3,5,0)
arrangeModes = c(TRUE, TRUE, FALSE)
continuousModes = c(FALSE,FALSE,TRUE)
parafac4microbiome::plotPARAFACmodel(cp_tongue$Fac, processedTongue, 2, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Tongue")
parafac4microbiome::plotPARAFACmodel(cp_lowling$Fac, processedLowling, 1, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Lowling")
parafac4microbiome::plotPARAFACmodel(cp_lowinter$Fac, processedLowinter, 2, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Lowinter")
parafac4microbiome::plotPARAFACmodel(cp_upling$Fac, processedUpling, 2, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Upling")
parafac4microbiome::plotPARAFACmodel(cp_upinter$Fac, processedUpinter, 2, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Upinter")
parafac4microbiome::plotPARAFACmodel(cp_saliva$Fac, processedSaliva, 2, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Saliva")
parafac4microbiome::plotPARAFACmodel(cp_metabolomics$Fac, processedMetabolomics, 1, colourCols, legendTitles, xLabels, legendColNums, arrangeModes, continuousModes,
overallTitle = "Metabolomics")
