问题描述：
因为我的数据每个时间点有多个重复样本，所以我想要使用circacompare包中的circacompare_mixed()比较两组（对照组和实验组）基因表达量的节律特征，并且基因数量达5000多个（这里只摘取了其中两个），因此我又用了for循环，但跑不出结果且出现了报错和警示：Iteration 1, LME step: nlminb() did not converge (code = 1). Do increase 'msMaxIter'!Error in nlme.formula(model = form_group, random = randomeffects_formula, :
Unable to form Cholesky decomposition: the leading minor of order 1 is not pos.def.
我的数据：
dput(t)
structure(list(id = c(1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3,
3, 3, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3), group = c("CON",
"CON", "CON", "CON", "CON", "CON", "CON", "CON", "CON", "CON",
"CON", "CON", "CON", "CON", "CON", "CC", "CC", "CC", "CC", "CC",
"CC", "CC", "CC", "CC", "CC", "CC", "CC", "CC", "CC", "CC"),
time = c(0, 6, 12, 18, 24, 0, 6, 12, 18, 24, 0, 6, 12, 18,
24, 0, 6, 12, 18, 24, 0, 6, 12, 18, 24, 0, 6, 12, 18, 24),
Gnai3 = c(17.51364, 17.0441, 20.73355, 16.02792, 9.2444,
20.62992, 17.28127, 29.32941, 14.5062, 16.66473, 8.196866,
8.362214, 39.94942, 16.29845, 19.80216, 13.74092, 14.08834,
14.04772, 26.00864, 17.1612, 14.19446, 12.62367, 13.53981,
43.95303, 18.80314, 15.35234, 16.63615, 23.49265, 39.86695,
16.80003), Scmh1 = c(1.151244, 1.007064, 0.8978803, 0.8906483,
1.073716, 1.133559, 0.9923863, 0.3706273, 0.9956849, 1.540639,
0.9482186, 1.105261, 0.1125411, 1.31828, 1.304933, 1.016301,
1.008043, 1.098772, 0.3921483, 1.404781, 1.214758, 1.258968,
0.8213903, 0.2767058, 1.430973, 0.9601138, 1.321019, 0.6648583,
0.3420466, 1.229267)), class = "data.frame", row.names = c(NA,
30L))
我的代码：
{ampm<-data.frame()
amp_1m<-data.frame()
mesor_1m<-data.frame()
phase_1m<-data.frame()
amp_2m<-data.frame()
mesor_2m<-data.frame()
phase_2m<-data.frame()
amp_pm<-data.frame()
mesor_pm<-data.frame()
phase_pm<-data.frame()
rhy_p1m<-data.frame()
rhy_p2m<-data.frame()
}
for(i in 4:5){colnames(t[i])->gm
out_im<-circacompare_mixed(x=t,col_time="time", col_outcome=gm,col_group ="group",col_id="id",control=list(grouped_params=c("phi","alpha","k"),random_params=c("phi1","alpha1","k1")),period = 24,alpha_threshold = 0.05)
out_im[[2]]$gene_symbol<-gm
parameter_im<-as.data.frame(out_im[[2]])
mesor1m<-slice(parameter_im,3)
mesor2m<-slice(parameter_im,4)
mesorpm<-slice(parameter_im,6)
amp1m<-slice(parameter_im,7)
amp2m<-slice(parameter_im,8)
amppm<-slice(parameter_im,10)
phase1m<-slice(parameter_im,11)
phase2m<-slice(parameter_im,12)
phasepm<-slice(parameter_im,14)
amp_1m<-rbind(amp_1m,amp1m)
mesor_1m<-rbind(mesor_1m,mesor1m)
phase_1m<-rbind(phase_1m,phase1m)
amp_2m<-rbind(amp_2m,amp2m)
mesor_2m<-rbind(mesor_2m,mesor2m)
phase_2m<-rbind(phase_2m,phase2m)
amp_pm<-rbind(amp_pm,amppm)
mesor_pm<-rbind(mesor_pm,mesorpm)
phase_pm<-rbind(phase_pm,phasepm)
}
我的运行结果：
Warning: Iteration 1, LME step: nlminb() did not converge (code = 1). Do increase 'msMaxIter'! PORT message: false convergence (8)Warning: Singular precision matrix in level -1, block 1Warning: Singular precision matrix in level -1, block . Error in nlme.formula(model = form_group, random = randomeffects_formula, :
Unable to form Cholesky decomposition: the leading minor of order 1 is not pos.def.
我的系统环境：
R version 4.2.1 (2022-06-23 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19044)
Matrix products: default
locale:
[1] LC_COLLATE=Chinese (Simplified)_China.utf8
[2] LC_CTYPE=Chinese (Simplified)_China.utf8
[3] LC_MONETARY=Chinese (Simplified)_China.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=Chinese (Simplified)_China.utf8
attached base packages:
[1] stats graphics grDevices utils datasets methods
[7] base
loaded via a namespace (and not attached):
[1] fgsea_1.24.0 colorspace_2.0-3
[3] ggtree_3.6.2 rjson_0.2.21
[5] gson_0.0.9 circlize_0.4.15
[7] qvalue_2.30.0 XVector_0.38.0
[9] GlobalOptions_0.1.2 aplot_0.1.9
[11] clue_0.3-64 rstudioapi_0.14
[13] farver_2.1.1 graphlayouts_0.8.4
[15] ggrepel_0.9.2 bit64_4.0.5
[17] AnnotationDbi_1.60.0 fansi_1.0.3
[19] scatterpie_0.1.8 codetools_0.2-19
[21] splines_4.2.1 doParallel_1.0.17
[23] cachem_1.0.6 GOSemSim_2.24.0
[25] knitr_1.42 polyclip_1.10-4
[27] jsonlite_1.8.4 cluster_2.1.4
[29] GO.db_3.16.0 png_0.1-8
[31] ggforce_0.4.1 compiler_4.2.1
[33] httr_1.4.5 assertthat_0.2.1
[35] Matrix_1.5-3 fastmap_1.1.0
[37] lazyeval_0.2.2 cli_3.6.0
[39] tweenr_2.0.2 htmltools_0.5.4
[41] tools_4.2.1 igraph_1.3.5
[43] gtable_0.3.1 glue_1.6.2
[45] GenomeInfoDbData_1.2.9 reshape2_1.4.4
[47] dplyr_1.0.10 fastmatch_1.1-3
[49] Rcpp_1.0.9 enrichplot_1.18.3
[51] Biobase_2.58.0 vctrs_0.5.1
[53] Biostrings_2.66.0 ape_5.6-2
[55] nlme_3.1-161 iterators_1.0.14
[57] ggraph_2.1.0 xfun_0.36
[59] stringr_1.5.0 lifecycle_1.0.3
[61] clusterProfiler_4.6.0 DOSE_3.24.2
[63] zlibbioc_1.44.0 MASS_7.3-58.1
[65] scales_1.2.1 tidygraph_1.2.2
[67] parallel_4.2.1 RColorBrewer_1.1-3
[69] yaml_2.3.7 ComplexHeatmap_2.14.0
[71] memoise_2.0.1 gridExtra_2.3
[73] ggplot2_3.4.1 downloader_0.4
[75] ggfun_0.0.9 HDO.db_0.99.1
[77] yulab.utils_0.0.6 stringi_1.7.12
[79] RSQLite_2.2.20 circacompare_0.1.1
[81] S4Vectors_0.36.1 foreach_1.5.2
[83] tidytree_0.4.2 BiocGenerics_0.44.0
[85] BiocParallel_1.32.5 shape_1.4.6
[87] GenomeInfoDb_1.34.9 rlang_1.0.6
[89] pkgconfig_2.0.3 matrixStats_0.63.0
[91] bitops_1.0-7 evaluate_0.20
[93] lattice_0.20-45 purrr_1.0.1
[95] treeio_1.22.0 patchwork_1.1.2
[97] cowplot_1.1.1 shadowtext_0.1.2
[99] bit_4.0.5 tidyselect_1.2.0
[101] plyr_1.8.8 magrittr_2.0.3
[103] R6_2.5.1 IRanges_2.32.0
[105] generics_0.1.3 DBI_1.1.3
[107] pillar_1.8.1 withr_2.5.0
[109] KEGGREST_1.38.0 RCurl_1.98-1.9
[111] tibble_3.1.8 crayon_1.5.2
[113] utf8_1.2.2 rmarkdown_2.20
[115] viridis_0.6.2 GetoptLong_1.0.5
[117] grid_4.2.1 data.table_1.14.6
[119] blob_1.2.3 digest_0.6.31
[121] tidyr_1.2.1 gridGraphics_0.5-1
[123] stats4_4.2.1 munsell_0.5.0
[125] viridisLite_0.4.1 ggplotify_0.1.0
