sceList=lapply(samples,function(pro){ #pro=samples[1] print(pro) tmp=Read10X(file.path(dir,pro)) if(length(tmp)==2){ ct=tmp[[1]] }else{ct=tmp} sce=CreateSeuratObject(counts=ct, project=pro, min.cells=5, min.features=300) return(sce) })

时间: 2025-03-23 07:15:19 浏览: 52
这段代码的主要功能是对一组样本 (`samples`) 进行循环处理,并对每个样本创建一个 Seurat 对象 (`sce`),最终将所有结果存储在一个列表 `sceList` 中。以下是逐步解析: --- ### **代码详解** 1. **函数框架** 使用了 R 语言中的 `lapply()` 函数,它会对输入向量 `samples` 的每一个元素应用指定的匿名函数。 2. **打印当前样本名** ```R print(pro) ``` 每次迭代时会打印出正在处理的样本名称 `pro`,方便调试或跟踪进度。 3. **读取数据文件** ```R tmp = Read10X(file.path(dir, pro)) ``` - 调用了 `Read10X()` 函数从指定路径加载数据。这里假设 `dir` 是存放数据的目录。 - `file.path(dir, pro)` 构造了一个完整的文件路径,指向某个特定样本的数据文件。 4. **检查数据结构并提取矩阵** ```R if (length(tmp) == 2) { ct = tmp[[1]] } else { ct = tmp } ``` - 如果返回值是一个长度为 2 的列表,则取出第一个元素作为计数矩阵 `ct`;否则直接将整个 `tmp` 视作计数矩阵。 5. **创建 Seurat 对象** ```R sce = CreateSeuratObject( counts = ct, project = pro, min.cells = 5, min.features = 300 ) ``` - 利用 Seurat 包提供的 `CreateSeuratObject()` 创建单细胞数据分析的对象。 - 参数说明: - `counts`: 计数矩阵; - `project`: 样本标识符(即 `pro`); - `min.cells`: 筛选掉表达在少于 5 个细胞中的基因; - `min.features`: 筛选掉含有少于 300 个特征的细胞。 6. **返回结果** 每次迭代都会返回生成的 Seurat 对象 `sce`,最后由 `lapply()` 将所有的结果收集到列表 `sceList` 中。 --- ### **用途及意义** 此段代码适用于批量处理单细胞测序数据集的情景,通过自动化的方式逐一读取样本、构建分析所需的 Seurat 对象集合。这对于后续统一进行质控过滤、归一化、降维聚类等步骤非常有用。 --- ###
阅读全文

相关推荐

library(nloptr) library(cluster) library(data.table) library(purrr) library(iterators) library(parallel) library(doParallel) registerDoParallel(cores = 7) # 根据实际CPU核心数调整 T <- 507 # 时间期数 K <- 94 # 特征数 G <- 6 # 组数 #=============================== # 0. 数据预处理(提前分割数据) #=============================== # 转换为数据框 df_standard <- as.data.frame(df_standard) # 假设 df_standard 包含 year_m, chara_name, ia_col, ret 等列 setDT(df_standard) # 转换为data.table加速 # 按 sample_ym 的顺序提取子数据表 split_data <- lapply(sample_ym, function(year) { df_standard[year_m == year] }) # 预计算每个时间点的公司数量N N <- df_standard[, .(N = .N), by = year_m][order(year_m)] #=============================== # 1. 生成数据列表(避免循环内重复筛选) #=============================== data_list <- map(split_data, ~{ X_t <- as.matrix(.x[, c(chara_name, ia_col), with = FALSE]) R_t <- .x$ret list(X = X_t, R = R_t) }) #=============================== # 2. 计算availability矩阵 #=============================== # 标记特征k在时间t是否可用(至少有一个非缺失值) availability <- matrix(FALSE, nrow = K, ncol = T) for (t in 1:T) { X_t <- data_list[[t]]$X availability[, t] <- apply(X_t, 2, function(col) any(!is.na(col))) } #=============================== # 3. 构造稀疏设计矩阵(高效实现) #=============================== # 预先生成所有块矩阵(使用replace快速处理NA) X_blocks <- map2(data_list, 1:T, function(data, t) { X_t <- data$X X_t[, !availability[, t]] <- 0 # 不可用特征列置0 X_t[is.na(X_t)] <- 0 # 缺失值置0 X_t }) # 直接构造块对角稀疏矩阵 X_block <- bdiag(X_blocks) X_block <- as(X_block, "CsparseMatrix") # 合并收益率向量 R_all <- unlist(map(data_list, "R")) #=============================== # 3. 初始化参数(处理缺失值) #=============================== lambda <- 5 gamma = t(gamma) #排除掉nan的值 gamma[is.na(gamma)] <- 0 alpha_g = t(alpha_g) #=============================== # 4. 定义目标函数与梯度(严格乘积项) #=============================== objective <- function(gamma_vec, alpha, availability, X_block, R_all, lambda) { gamma <- matrix(gamma_vec, nrow = K, ncol = T) pred <- as.vector(X_block %*% c(gamma)) # 自动跳过缺失值(对应位置为0) error <- sum((R_all - pred)^2) penalty <- 0 for (k in 1:1) { valid_t <- which(availability[k, ]) if (length(valid_t) == 0) next # 跳过全缺失特征 product_term <- 1 for (g in 1:1) { gamma_k_valid <- gamma[k, valid_t] alpha_g_valid <- alpha[g, valid_t] product_term <- product_term * sqrt(sum((gamma_k_valid - alpha_g_valid)^2 + 1e-8)) } penalty <- penalty + product_term } obj = error + lambda * penalty print(obj) return(obj) } gradient <- function(gamma_vec, alpha, availability, X_block, R_all, lambda) { print('ok') K <- nrow(availability) T <- ncol(availability) G <- nrow(alpha) gamma <- matrix(gamma_vec, nrow = K, ncol = T) # 预计算预测误差项的梯度 pred <- as.vector(X_block %*% gamma_vec) error_grad <- (2 / length(R_all)) * crossprod(X_block, (pred - R_all)) # 预计算惩罚项的梯度 penalty_grad <- matrix(0, nrow = K, ncol = T) for (k in 1:K) { valid_t <- which(availability[k, ]) if (length(valid_t) == 0) next gamma_k <- gamma[k, valid_t] norms <- sapply(1:G, function(g) { sqrt(sum((gamma_k - alpha[g, valid_t])^2) + 1e-8) }) product <- prod(norms) for (g in 1:G) { term <- (gamma_k - alpha[g, valid_t]) / norms[g] term <- term * (product / norms[g]) penalty_grad[k, valid_t] <- penalty_grad[k, valid_t] + term } } penalty_grad <- lambda * penalty_grad # 合并梯度 total_grad <- error_grad + c(penalty_grad) return(total_grad) } #=============================== # 5. 交替优化主循环 #=============================== max_iter <- 20 tol <- 1e-1 converged <- FALSE for (iter in 1:max_iter) { #记录次数 print(iter) gamma_old <- gamma alpha_old <- alpha_g # --- 步骤1:固定α,优化γ --- opts <- list( algorithm = "NLOPT_LD_MMA", # 更适合非凸问题的算法 xtol_rel = 1e-4, # 降低收敛精度要求 maxeval = 10 # 减少每次迭代的最大评估次数 ) result <- nloptr( x0 = c(gamma), eval_f = function(x) objective(x, alpha_g, availability, X_block, R_all, lambda), #梯度的方向 eval_grad_f = function(x) gradient(x, alpha_g, availability, X_block, R_all, lambda), opts = opts ) gamma <- matrix(result$solution, nrow = K, ncol = T) # --- 步骤2:固定γ,更新α和分组 --- dist_matrix <- foreach(k = 1:K, .combine = "rbind") %dopar% { valid_t <- which(availability[k, ]) if (length(valid_t) == 0) return(rep(Inf, G)) gamma_k_valid <- gamma[k, valid_t] sapply(1:G, function(g) { alpha_g_valid <- alpha[g, valid_t] sqrt(sum((gamma_k_valid - alpha_g_valid)^2)) }) } groups <- apply(dist_matrix, 1, which.min) # alpha更新 alpha_g <- matrix(NA, nrow = G, ncol = T) for (g in 1:G) { members <- which(groups == g) if (length(members) > 0) { alpha_g[g, ] <- colMeans(gamma[members, , drop = FALSE], na.rm = TRUE) } } delta_gamma <- max(abs(gamma - gamma_old), na.rm = TRUE) delta_alpha <- max(abs(alpha_g - alpha_old), na.rm = TRUE) cat(sprintf("Iter %d: Δγ=%.4f, Δα=%.4f\n", iter, delta_gamma, delta_alpha)) if (delta_gamma < tol && delta_alpha < tol) break }Error in nloptr(x0 = c(gamma), eval_f = function(x) objective(x, alpha_g, : REAL() can only be applied to a 'numeric', not a 'S4'

# 随机森林多分类模型 # packages 安装 —— install.packages("") setwd("D:/Rpackages") # packages 载入 —— library() library(randomForest) library(tidyverse) library(skimr) library(DataExplorer) library(caret) library(pROC) library(ggplot2) library(splines) library(reshape2) library(scales) library(ggprism) library(ggpubr) # 加载数据 boston <- read.csv(file.choose()) # 数据框转换 boston <- as.data.frame(boston) # boston <- na.pass(boston) # boston <- na.omit(boston) # na.action = na.omit 返回删除所有包含NA值的行的数据框或矩阵 # na.action = na.pass 将NA值传递给后续的计算,可能会导致结果中也出现NA # 查看数据概貌 skim(boston) # 查看数据缺失情况 plot_missing(boston) boston[is.na(boston)] <- 0 #将matrix(此处为boston)中的NA值替换为0,按需运行 # 数据类型修正,将boston中的i列转为factor,i 将依次取c()中的值 for (i in c(1:2)) { boston[,i] <- factor(boston[,i]) } # 因变量分布情况 table(boston$Class) # 拆分数据 set.seed(1) trains <- createDataPartition( y = boston$Class, p = 0.7, list = F ) traindata <- boston[trains,] testdata <- boston[-trains,] # # test:自选训练集测试集(无特殊可不理 # traindata <- read.csv(file.choose()) # traindata <- as.data.frame(traindata) # for (i in c(1:5)) { # traindata[,i] <- factor(traindata[,i]) # } # # testdata <- read.csv(file.choose()) # testdata <- as.data.frame(testdata) # for (i in c(1:5)) { # testdata[,i] <- factor(testdata[,i]) # } # 拆分后因变量分布 table(traindata$Class) table(testdata$Class) dim(boston) # # 构建公式,自变量列名~a到b列(colnames[a:b]) # # colnames(boston) 用于获取数据框的列名, # form_clsm <- as.formula( # paste0( # "class~", # paste(colnames(traindata)[3:112],collapse = "+") # ) # ) # form_clsm # # # 构建模型 # set.seed(100001) # fit_rf_clsm <- randomForest( # form_clsm, # data = traindata, # ntree = 500, # 默认500,error & trees不稳定需增加树的数量 # mtry = 6, # 每个节点可提供选择的变量数目(指定节点中用于二叉树的变量个数 # # 默认情况下数据集变量个数的二次方根(分类模型)或三分之一(预测模型)) # importance = T # importance = TRUE: 启用特征重要性计算 # # 随机森林会根据每个特征对模型性能的影响,计算出每个特征的重要性 # ) # print(fit_rf_clsm) # 定义训练集特征和目标变量 X_train <- traindata[, -(1:2)] #traindata中除第1、2列外的列为自变量 x y_train <- as.factor(traindata[, 2]) #traindata中第2列为因变量 y # # 创建随机森林分类模型(基准模型 # model <- randomForest(x = X_train, y = y_train, ntree = 500) # # # 创建训练控制对象 # ctrl <- trainControl(method = "repeatedcv", number = 5, repeats = 10) # # k 折交叉验证作为模型评估方法,小数据做重复10次5折交叉验证 # # # 进行参数调优 # # mtry参数调节 # grid1 <- expand.grid(mtry = c(48:52)) # 定义mtry范围 # # grid1 <- expand.grid(mtry = c(10, 20, 30, 40, 50, 60, 70, 80, 90, 100)) # # 进阶 # # num_features <- ncol(X_train) # # mtry_range <- floor(sqrt(num_features)) # 取平方根并向下取整 # # grid_mtry <- expand.grid(mtry = seq(max(2, mtry_range - 5), # # min(mtry_range + 5, num_features))) # # # 使用caret包进行调参 # rf_model <- train(x = X_train, y = y_train, # 自变量 x, 因变量 y # method = "rf", # 使用的模型方法为随机森林(rf) # trControl = ctrl, # tuneGrid = grid1,) # # 输出结果 # print(rf_model) # # # # 定义最佳mtry参数 # grid2 <- expand.grid(mtry = c(30)) # # # 定义模型列表,存储每一个模型评估结果 # modellist <- list() # # 调整的参数是决策树的数量 # for (ntree in seq(100, 2000, by=100)) { # seq函数构建一个集合,间距by = 100 # set.seed(101) # fit <- train(x = X_train, y = y_train, method="rf", # metric="Accuracy", tuneGrid=grid2, # trControl=ctrl, ntree=ntree) # key <- toString(ntree) # modellist[[key]] <- fit # } # # # compare results # results <- resamples(modellist) # # # 输出最佳模型和参数 # summary(results) # # # 可视化模型性能,accuracy越近1越好,kappa越近1越好 # bwplot(results) # 箱线图 # dotplot(results) # 点图 # densityplot(results) # 密度图 # 使用最佳参数训练模型 set.seed(1001) fit_rf_clsm <- randomForest(x = X_train, y = y_train, mtry = 50, ntree = 500, importance = T) print(fit_rf_clsm) # ntree参数与error之间的关系图示 plot(fit_rf_clsm,main = "ERROR & TREES") legend("top", legend = colnames(fit_rf_clsm$err.rate), lty = 1:6, col = 1:6, horiz = T, cex = 0.9) plot(randomForest::margin(fit_rf_clsm), main = '观测值被判断正确的概率图') # 变量重要性 varImpPlot(fit_rf_clsm,main ="varImpPlot") varImpPlot(fit_rf_clsm,main = "varImpPlot",type = 1) varImpPlot(fit_rf_clsm,main = "varImpPlot",type = 2) importance_genus <- data.frame(importance(fit_rf_clsm)) importance_genus <- importance_genus[order(importance_genus$MeanDecreaseGini, decreasing=TRUE),] importance_genus <- importance_genus[order(importance_genus$MeanDecreaseAccuracy, decreasing=TRUE),] head(importance_genus) write.table(importance_genus,"importance_genus_DaChuang_SiteClass_unexposed20250209.txt", sep = '\t',col.names = NA,quote = FALSE) # 保存重要特征到文件 # # 预测 # # 训练集预测概率 # trainpredprob <- predict(fit_rf_clsm,newdata = traindata,type = "prob") # # 训练集ROC # multiclass.roc(response = traindata$class,predictor = trainpredprob) # # 训练集预测分类 # trainpredlab <- predict(fit_rf_clsm,newdata = traindata,type = "Class") # # 训练集混淆矩阵 # confusionMatrix_train <- confusionMatrix(data = trainpredlab, # reference = traindata$Class, # mode = "everything") # # 训练集综合结果 # multiClassSummary( # data.frame(obs = traindata$Class,pred=trainpredlab), # lev = levels(traindata$Class) # ) # # # # 作图展示 top N 重要的 OTUs # varImpPlot(fit_rf_clsm, n.var = min(20, nrow(fit_rf_clsm$importance)), # main = 'Top 20 - variable importance') # # # 测试集预测概率 # # X_test <- testdata[, -(1:2)] #testdata中除第1、2列外的列为自变量 x # # y_test <- as.factor(testdata[, 2]) #testdata中第2列为因变量 y # testpredprob <- predict(fit_rf_clsm,newdata = testdata,type = "prob") # # 测试集ROC # multiclass.roc(response = testdata$Class,predictor = testpredprob) # # 测试集预测分类 # testpredlab <- predict(fit_rf_clsm,newdata = testdata,type = "Class") # # 测试集混淆矩阵 # confusionMatrix_test <- confusionMatrix(data = testpredlab, # reference = testdata$Class, # mode = "everything") # confusionMatrix_test # # 测试集综合结果 # multiClassSummary( # data.frame(obs = testdata$Class,pred=testpredlab), # lev = levels(testdata$Class) # ) # 交叉验证帮助选择特定数量的特征 # 5次重复十折交叉验证 set.seed(10001) otu_train.cv <- replicate(5, rfcv(traindata[,-(1:2)], # 除去第a到b列(a:b) traindata$Class, cv.fold = 10, step = 1.5), simplify = FALSE) otu_train.cv <- data.frame(sapply(otu_train.cv, '[[', 'error.cv')) otu_train.cv$otus <- rownames(otu_train.cv) otu_train.cv <- reshape2::melt(otu_train.cv, id = 'otus') otu_train.cv$otus <- as.numeric(as.character(otu_train.cv$otus)) otu_train.cv.mean <- aggregate(otu_train.cv$value, by = list(otu_train.cv$otus), FUN = mean) head(otu_train.cv.mean, 18) # 绘图观察拐点 p <- ggplot(otu_train.cv,aes(otus,value)) + geom_smooth(se = FALSE, method = 'glm',formula = y~ns(x,6)) + theme(panel.grid = element_blank(), panel.background = element_rect(color = 'black',fill = 'transparent')) + labs(title = '',x='Number of genus',y='Cross-validation error') p # 在横坐标(xintecept)绘制竖线 p + geom_vline(xintercept = 160) # # 备用 # p2 <- ggplot(otu_train.cv,aes(otus,value)) + # geom_line() + # theme(panel.grid = element_blank(), # panel.background = element_rect(color = 'black', fill = 'transparent')) + # labs(title = '',x = 'Number of OTUs', y = 'Cross-validation error') # p2 # # p2 + geom_vline(xintercept = 30) # 大约提取前 N个重要的特征 importance_genus[1:160,] # importance_Ngenus <- importance_genus[1:160,] # # 输出表格 # write.table(importance_genus[1:70, ], # 'importance_genus_top70_of_zhangxiaofeng_human_drowning.txt', # sep = '\t', col.names = NA, quote = FALSE) # # 变量重要性 # varImpPlot(fit_rf_clsm,main ="varImpPlot") # varImpPlot(fit_rf_clsm,main = "varImpPlot",type = 1) # varImpPlot(fit_rf_clsm,main = "varImpPlot",type = 2) # # varImpPlot(fit_rf_clsm, n.var = min(160, nrow(fit_rf_clsm$importance)), # min(num,)处的num为图中的菌属数量 # main = 'Top 19 - variable importance',type = 1) # 简约分类器(只取部分高影响的自变量) # 选择 top N 重要的 OTUs,例如上述已经根据“Mean Decrease Accuracy”排名获得 genus_select <- rownames(importance_genus)[1:160] # 数据子集的训练集和测试集 genus_train_top <- traindata[ ,c(genus_select, 'Class')] genus_test_top<- testdata[ ,c(genus_select, 'Class')] # set.seed(10001) # form_clsm1 <- as.formula( # paste0( # "class~", # paste(colnames(genus_train_top)[1:10],collapse = "+") # ) # ) # 构建模型 # fit_rf_clsm1 <- randomForest( # form_clsm1, # data = genus_train_top, # ntree = 500, # mtry = 6, # importance = T # ) x_train1 <- genus_train_top[, -160 - 1] # 自变量 x y_train1 <- as.factor(genus_train_top[, 160 + 1]) # 因变量 y # fit_rf_clsm_test1 <- randomForest(x = x_train1, # y = y_train1, # ntree = 500, # 增加树的数量以提高稳定性 # importance = TRUE # 启用特征重要性计算 # ) # # fit_rf_clsm_test1 # # # 5 折交叉验证,重复 10 次 # ctrl1 <- trainControl(method = "repeatedcv", number = 5, repeats = 10) # # # 定义 mtry 和 ntree 的参数范围 # grid3 <- expand.grid(mtry = c(2:15)) # 定义mtry范围 # # # 进阶 # # num_features <- ncol(X_train) # # mtry_range <- floor(sqrt(num_features)) # 取平方根并向下取整 # # grid_mtry <- expand.grid(mtry = seq(max(2, mtry_range - 5), # # min(mtry_range + 5, num_features))) # # # 使用caret包进行调参 # fit_rf_clsm_test2 <- train(x = x_train1, y = y_train1, # 自变量 x, 因变量 y # method = "rf", # 使用的模型方法为随机森林(rf) # trControl = ctrl1, # tuneGrid = grid3,) # # 输出结果 # print(fit_rf_clsm_test2) # # # 定义最佳mtry参数 # grid4 <- expand.grid(mtry = c(16)) # # # 定义模型列表,存储每一个模型评估结果 # modellist1 <- list() # # 调整的参数是决策树的数量 # for (ntree in seq(100, 2000, by=100)) { # seq函数构建一个集合,间距by = 100 # set.seed(100003) # fit1 <- train(x = x_train1, y = y_train1, method="rf", # metric="Accuracy", tuneGrid=grid4, # trControl=ctrl1, ntree=ntree) # key1 <- toString(ntree) # modellist1[[key1]] <- fit1 # } # # # compare results # results1 <- resamples(modellist1) # # # 输出最佳模型和参数 # summary(results1) # # # 可视化模型性能,accuracy越近进1越好,kappa越近1越好 # bwplot(results1) # 箱线图 # dotplot(results1) # 点图 # densityplot(results1) # 密度图 # 使用最佳参数训练模型 set.seed(1) fit_rf_clsm1 <- randomForest(x = x_train1, y = y_train1, mtry = 12, ntree = 500, importance = T) print(fit_rf_clsm1) # # "ERROR & TREES" # plot(fit_rf_clsm1,main = "ERROR & TREES") # # # plot(randomForest::margin(fit_rf_clsm1), main = '观测值被判断正确的概率图') # # 预测 # # 训练集预测概率 # trainpredprob <- predict(fit_rf_clsm1,newdata = genus_train_top,type = "prob") # # 训练集ROC # multiclass.roc(response = genus_train_top$Class,predictor = trainpredprob) # # 训练集预测分类 #trainpredlab <- predict(fit_rf_clsm1,newdata = genus_train_top,type = "Class") # # 训练集混淆矩阵 # confusionMatrix(data = trainpredlab, # reference = genus_train_top$Class, # mode = "everything") # 预测 # 测试集预测概率 # x_test1 <- genus_test_top[, -11] # 自变量 x # y_test1 <- as.factor(genus_test_top[, 11]) # 因变量 y testpredprob <- predict(fit_rf_clsm1,newdata = genus_test_top,type = "prob") write.table(testpredprob, file = "D:/Rpackages/testpredprob.txt", sep = "\t", row.names = FALSE, col.names = TRUE) # 测试集ROC multiclass.roc(response = genus_test_top$Class,predictor = testpredprob) # 测试集预测分类 testpredlab <- predict(fit_rf_clsm1,newdata = genus_test_top,type = "Class") # 测试集混淆矩阵 confusion_matrix <- confusionMatrix(data = testpredlab, reference = genus_test_top$Class, mode = "everything") # 测试集综合结果 multiClassSummary( data.frame(obs = genus_test_top$Class,pred=testpredlab), lev = levels(genus_test_top$Class) ) # 查看样本预测结果 results <- data.frame(Actual = genus_test_top$Class, Predicted = testpredlab) # 测试集预测分类 # testpredlab <- predict(fit_rf_clsm1,newdata = testdata,type = "class") # t <- table(testpredlab,testdata$class) # acc = sum(diag(t))/nrow(testdata)*100 # print(paste("模型准确率为:",round(acc,4),sep='')) # 绘制混淆矩阵热图(内容复杂,好汉谨慎处之) # confusion_matrix是混淆矩阵对象 # 转换混淆矩阵为数据框 roc1 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,1:6]) roc1 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,1]) plot(roc1$rocs[[1]],col="#1f77b4",print.auc = TRUE,print.auc.x=0.8,print.auc.y=0.8) roc2 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,2]) plot(roc2$rocs[[1]],add=TRUE,col="#ff7f0e",print.auc = TRUE,print.auc.x=0.6,print.auc.y=0.6) roc3 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,3]) plot(roc3$rocs[[1]],add=TRUE,col="#2ca02c",print.auc=TRUE,print.auc.x=0.5,print.auc.y=0.5) roc4 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,4]) plot(roc4$rocs[[1]],add=TRUE,col="#d62728",print.auc=TRUE,print.auc.x=0.4,print.auc.y=0.4) roc5 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,5]) plot(roc5$rocs[[1]],add=TRUE,col="#9467bd",print.auc=TRUE,print.auc.x=0.3,print.auc.y=0.3) roc6 <- multiclass.roc(response = genus_test_top$Class,predictor =testpredprob[,6]) plot(roc1$rocs[[6]], add = TRUE, col = "#8c564b", print.auc = TRUE, print.auc.x = 0.2, print.auc.y = 0.2) # confusion_matrix_df <- as.data.frame.matrix(confusion_matrix$table) # colnames(confusion_matrix_df) <- c("F","H") # rownames(confusion_matrix_df) <- c("F","H") # # #c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15") # colnames(confusion_matrix_df) <- c("01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12", "13", "14", "15") # rownames(confusion_matrix_df) <- c("01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12", "13", "14", "15") # # 计算归一化值 # draw_data <- round(confusion_matrix_df / colSums(confusion_matrix_df),2) # draw_data <- as.data.frame(t(draw_data)) #使用t()函数转换draw_data的行列,原本列为真实类别,现转为行 # draw_data$real <- rownames(draw_data) #提取行名作为新生成列real的类别 # draw_data <- melt(draw_data) #将宽格式的数据框转换为长格式 # # # 绘制矩阵热图 # confusion1 <- ggplot(draw_data, aes(real, variable, fill = value)) + # geom_tile() + # geom_text(aes(label = scales::percent(value))) + # scale_fill_gradient(low = "#F0F0F0", high = "#3575b5") + # labs(x = "Prediction", y = "Reference", title = "Confusion matrix") + # theme_prism(border = T) + # theme(panel.border = element_blank(), # axis.ticks.y = element_blank(), # axis.ticks.x = element_blank(), # legend.position="right", # plot.title = element_text(hjust = 0.5)) # # confusion1 # # 二分类ROC曲线绘制 # # 计算ROC曲线的参数 # roc_obj <- roc(response = genus_test_top$class, predictor = testpredprob[, 2]) # roc_auc <- auc(roc_obj) # # 将ROC对象转换为数据框 # roc_data <- data.frame(1 - roc_obj$specificities, roc_obj$sensitivities) # # 绘制ROC曲线 # ROC1 <- ggplot(roc_data, aes(x = 1 - roc_obj$specificities, y = roc_obj$sensitivities)) + # geom_line(color = "#0073C2FF", size = 1) + # annotate("segment", x = 0, y = 0, xend = 1, yend = 1, linetype = "dashed", color = "gray") + # annotate("text", x = 0.8, y = 0.2, label = paste("AUC =", round(roc_auc, 3)), size = 4, color = "black") + # coord_cartesian(xlim = c(0, 1), ylim = c(0, 1)) + # theme_pubr() + # labs(x = "1 - Specificity", y = "Sensitivity") + # ggtitle("ROC Curve") + # theme(plot.title = element_text(size = 14, face = "bold")) + # theme_prism(border = T) # ROC1 # # # 计算 ROC 和 AUC # roc_obj <- roc(response = genus_test_top$class, predictor = testpredprob[, 2]) # roc_auc <- auc(roc_obj) # # # 将 ROC 对象转换为数据框 # roc_data <- data.frame( # FPR = 1 - roc_obj$specificities, # TPR = roc_obj$sensitivities # ) # # # 平滑处理 # smooth_roc <- data.frame( # FPR = spline(roc_data$FPR, n = 500)$x, # TPR = spline(roc_data$TPR, n = 500)$y # ) # # # 绘制平滑后的 ROC 曲线 # ROC1 <- ggplot(smooth_roc, aes(x = FPR, y = TPR)) + # geom_line(color = "#0073C2FF", size = 1) + # annotate("segment", x = 0, y = 0, xend = 1, yend = 1, linetype = "dashed", color = "gray") + # annotate("text", x = 0.8, y = 0.2, label = paste("AUC =", round(roc_auc, 2)), size = 4, color = "black") + # coord_cartesian(xlim = c(0, 1), ylim = c(0, 1)) + # theme_pubr() + # labs(x = "1 - Specificity", y = "Sensitivity") + # ggtitle("Smoothed ROC Curve") + # theme(plot.title = element_text(size = 14, face = "bold")) + # theme_prism(border = T) # # ROC1 # # geom_smooth(se = FALSE, method = 'glm',formula = y~ns(x,6)) # # # 保存模型 # saveRDS(fit_rf_clsm, # file = "D:/Documents/R.data/fit_rf_clsm1_UnSimplifiedSiteClass_ConcernGender_UnexposedTop19_240209.rds") # # # 读取模型 # fit_rf_clsm1 <- readRDS("D:/Documents/R.data/fit_rf_clsm1_hand_Simplified_240102.rds") # # # 读取待分类数据 # testdata1 <- read.csv(file.choose()) # testdata1 <- as.data.frame(testdata) # for (i in c(1:2)) { # testdata1[,i] <- factor(testdata1[,i]) # } # # # 应用 # # 待分类数据集预测概率 # testpredprob <- predict(fit_rf_clsm1, newdata = testdata1, type = "prob") # # 测试集ROC # multiclass.roc(response = testdata1$class, predictor = testpredprob) # # 待分类数据集预测分类 # testpredlab <- predict(fit_rf_clsm1, newdata = testdata1,type = "class") # # 待分类数据集混淆矩阵 # confusion_matrix <- confusionMatrix(data = testpredlab, # reference = testdata1$class, # mode = "everything") # # 待分类数据集综合结果 # multiClassSummary( # data.frame(obs = testdata1$class,pred=testpredlab), # lev = levels(testdata1$class) # ) # # # 查看样本分类结果 # results <- data.frame(Actual = testdata1$class, Predicted = testpredlab) (这是整个的代码,testpredprob <- predict(fit_rf_clsm1,newdata = genus_test_top[1:5, ],type = "prob"),跑不下来)

Error in getGlobalsAndPackages(expr, envir = envir, globals = globals): The total size of the 10 globals exported for future expression (‘FUN()’) is 561.48 MiB.. This exceeds the maximum allowed size of 500.00 MiB (option 'future.globals.maxSize'). The three largest globals are ‘FUN’ (421.11 MiB of class ‘function’), ‘object’ (140.21 MiB of class ‘S4’) and ‘as.matrix’ (105.36 KiB of class ‘function’) Traceback: 1. NormalizeData.Seurat(st_data, normalization.method = "CLR", margin = 2, . notification = function(p) pb$update(p)) 2. NormalizeData(object = object[[assay]], normalization.method = normalization.method, . scale.factor = scale.factor, verbose = verbose, margin = margin, . ...) 3. NormalizeData.StdAssay(object = object[[assay]], normalization.method = normalization.method, . scale.factor = scale.factor, verbose = verbose, margin = margin, . ...) 4. NormalizeData(object = LayerData(object = object, layer = l, . fast = NA), normalization.method = normalization.method, . scale.factor = scale.factor, margin = margin, verbose = verbose, . ...) 5. NormalizeData.V3Matrix(object = LayerData(object = object, layer = l, . fast = NA), normalization.method = normalization.method, . scale.factor = scale.factor, margin = margin, verbose = verbose, . ...) 6. future_lapply(X = 1:ncol(x = chunk.points), FUN = function(i) { . block <- chunk.points[, i] . data <- if (margin == 1) { . object[block[1]:block[2], , drop = FALSE] . } . else { . object[, block[1]:block[2], drop = FALSE] . } . clr_function <- function(x) { . return(log1p(x = x/(exp(x = sum(log1p(x = x[x > 0]), . na.rm = TRUE)/length(x = x))))) . } . args <- list(data = data, scale.factor = scale.factor, verbose = FALSE, . custom_function = clr_function, margin = margin) . args <- args[names(x = formals(fun = norm.function))] . return(do.call(what = nor

大家在看

recommend-type

大唐杯仿真介绍.zip

大唐杯仿真 大唐杯仿真通常涉及通信网络的虚拟实践,特别是5G技术的相关应用。这类仿真旨在提供一个实践平台,让学生和参赛者能够在实际操作中深入理解和应用通信技术知识。 在大唐杯仿真中,参赛者可能会遇到多种任务和挑战,包括但不限于网络规划、设备配置、性能优化等。例如,在5G工程实践中,参赛者需要配置射频单元、光纤、光模块、电源类型等,这些都需要对5G通信技术有深入的了解。此外,车联网的仿真也是大唐杯的一个重点,参赛者需要配置车辆、路灯等模块,实现紧急前向防碰撞预警应用等功能。 大唐杯仿真通常在大赛平台(如学唐OnLine、虚拟仿真平台)上开放,供参赛者学习和训练。通过仿真实践,参赛者可以更加直观地了解通信网络的运行原理,提升实际操作能力,并锻炼解决实际问题的能力。
recommend-type

桌面便签_SimpleStickyNotes.zip

桌面便签_SimpleStickyNotes
recommend-type

美敦力BIS监护仪串口通讯协议手册

Document Title: BIS, MONITORING SYSTEMS, SERIAL PORT TECHNICAL SPEC
recommend-type

PFC与Fipy耦合技术:基于三角网格单元的双向流固耦合双轴压缩模拟,基于PFC流固耦合原理的双向耦合模拟技术:PFC与Fipy结合,三角网格单元实现渗流与双轴压缩模拟的双向交互作用 ,PFC流固耦合

PFC与Fipy耦合技术:基于三角网格单元的双向流固耦合双轴压缩模拟,基于PFC流固耦合原理的双向耦合模拟技术:PFC与Fipy结合,三角网格单元实现渗流与双轴压缩模拟的双向交互作用。,PFC流固耦合 PFC与Fipy结合,采用三角网格单元,双向耦合,实现渗流作用下的双轴压缩模拟。 ,PFC流固耦合; PFC与Fipy结合; 三角网格单元; 双向耦合; 渗流作用; 双轴压缩模拟。,PFC-Fipy流固双向耦合双轴压缩模拟
recommend-type

KR C4 电气元器件检修服务 - 系统技术.pdf

KR C4 电气元器件检修服务 - 系统技术.pdf

最新推荐

recommend-type

[学习]电子支付与网络银行第六章电子支付法律法规.ppt

[学习]电子支付与网络银行第六章电子支付法律法规.ppt
recommend-type

90%会计都不知道的Excel最常用九个快捷键!【2017-2018最新会计实务】.doc

90%会计都不知道的Excel最常用九个快捷键!【2017-2018最新会计实务】.doc
recommend-type

langchain4j-0.9.0.jar中文文档.zip

1、压缩文件中包含: 中文文档、jar包下载地址、Maven依赖、Gradle依赖、源代码下载地址。 2、使用方法: 解压最外层zip,再解压其中的zip包,双击 【index.html】 文件,即可用浏览器打开、进行查看。 3、特殊说明: (1)本文档为人性化翻译,精心制作,请放心使用; (2)只翻译了该翻译的内容,如:注释、说明、描述、用法讲解 等; (3)不该翻译的内容保持原样,如:类名、方法名、包名、类型、关键字、代码 等。 4、温馨提示: (1)为了防止解压后路径太长导致浏览器无法打开,推荐在解压时选择“解压到当前文件夹”(放心,自带文件夹,文件不会散落一地); (2)有时,一套Java组件会有多个jar,所以在下载前,请仔细阅读本篇描述,以确保这就是你需要的文件。 5、本文件关键字: jar中文文档.zip,java,jar包,Maven,第三方jar包,组件,开源组件,第三方组件,Gradle,中文API文档,手册,开发手册,使用手册,参考手册。
recommend-type

财务软件操作流程(1).doc

财务软件操作流程(1).doc
recommend-type

Visual-Studio-Team-System-XXXX-中的敏捷规划工具(1).docx

Visual-Studio-Team-System-XXXX-中的敏捷规划工具(1).docx
recommend-type

Wamp5: 一键配置ASP/PHP/HTML服务器工具

根据提供的文件信息,以下是关于标题、描述和文件列表中所涉及知识点的详细阐述。 ### 标题知识点 标题中提到的是"PHP集成版工具wamp5.rar",这里面包含了以下几个重要知识点: 1. **PHP**: PHP是一种广泛使用的开源服务器端脚本语言,主要用于网站开发。它可以嵌入到HTML中,从而让网页具有动态内容。PHP因其开源、跨平台、面向对象、安全性高等特点,成为最流行的网站开发语言之一。 2. **集成版工具**: 集成版工具通常指的是将多个功能组合在一起的软件包,目的是为了简化安装和配置流程。在PHP开发环境中,这样的集成工具通常包括了PHP解释器、Web服务器以及数据库管理系统等关键组件。 3. **Wamp5**: Wamp5是这类集成版工具的一种,它基于Windows操作系统。Wamp5的名称来源于它包含的主要组件的首字母缩写,即Windows、Apache、MySQL和PHP。这种工具允许开发者快速搭建本地Web开发环境,无需分别安装和配置各个组件。 4. **RAR压缩文件**: RAR是一种常见的文件压缩格式,它以较小的体积存储数据,便于传输和存储。RAR文件通常需要特定的解压缩软件进行解压缩操作。 ### 描述知识点 描述中提到了工具的一个重要功能:“可以自动配置asp/php/html等的服务器, 不用辛辛苦苦的为怎么配置服务器而烦恼”。这里面涵盖了以下知识点: 1. **自动配置**: 自动配置功能意味着该工具能够简化服务器的搭建过程,用户不需要手动进行繁琐的配置步骤,如修改配置文件、启动服务等。这是集成版工具的一项重要功能,极大地降低了初学者的技术门槛。 2. **ASP/PHP/HTML**: 这三种技术是Web开发中常用的组件。ASP (Active Server Pages) 是微软开发的服务器端脚本环境;HTML (HyperText Markup Language) 是用于创建网页的标准标记语言;PHP是服务器端脚本语言。在Wamp5这类集成环境中,可以很容易地对这些技术进行测试和开发,因为它们已经预配置在一起。 3. **服务器**: 在Web开发中,服务器是一个运行Web应用程序并响应客户端请求的软件或硬件系统。常见的服务器软件包括Apache、Nginx等。集成版工具提供了一个本地服务器环境,使得开发者可以在本地测试他们的应用程序。 ### 标签知识点 标签中仅出现了“PHP”一个关键词,这意味着该工具专注于与PHP相关的开发环境配置。 ### 压缩包子文件的文件名称列表知识点 1. **wamp.exe**: 这是Wamp5集成版工具的可执行文件,用户通过运行这个文件,即可启动Wamp5环境,开始进行PHP等相关开发。 2. **使用说明文档.txt**: 通常这样的文本文件包含了软件的安装指南、功能介绍、常见问题解答以及操作教程等内容。它是用户了解和使用软件的重要辅助资料,可以帮助用户快速上手。 ### 总结 Wamp5作为一个PHP集成版工具,为PHP开发者提供了一个便捷的本地开发环境。它通过自动化配置减少了开发者在搭建环境时的麻烦,整合了Web开发的关键组件,从而使得开发者能够专注于代码的编写,而不是环境配置。此外,它还提供了详尽的使用文档,帮助用户理解和运用这一工具,有效地提高了开发效率。对于希望快速入门PHP或需要一个简单易用的本地开发环境的开发者来说,Wamp5是一个非常好的选择。
recommend-type

【数据融合技术】:甘肃土壤类型空间分析中的专业性应用

# 摘要 数据融合技术作为一种集成多源数据信息的方法,在土壤类型空间分析中发挥着关键作用。本文介绍了数据融合技术的基本概念及其理论基础,阐述了数据预处理、同步整合及冲突解决等关键技术,并详细描述了甘肃土壤类型数据准备的流程,包括数据采集、质量评估、空间化处理及融合实践准备。通过具体案例分析,展示了数据融合在土壤类型空间分布分析、土壤质量评估及土壤保护规划中的应用。同时,文
recommend-type

sht20温湿度传感器使用什么将上拉电阻和滤波电容引出

<think>好的,我需要解决用户关于SHT20温湿度传感器的上拉电阻和滤波电容的引出方法问题。首先,我得回忆一下SHT20的基本特性和常见电路设计。SHT20是数字传感器,通常使用I2C接口,所以上拉电阻是必须的,用于确保数据线和时钟线的稳定。根据引用[2],SHT10是通过SCK和DATA线与单片机通信,而SHT30在引用[3]中使用I2C协议,需要上拉电阻。虽然用户问的是SHT20,但SHT系列通常设计类似,所以可以推断SHT20也需要类似的上拉电阻配置。通常I2C总线的上拉电阻值在4.7kΩ到10kΩ之间,但具体值可能取决于总线速度和电源电压。需要确认数据手册中的推荐值,但用户可能没有
recommend-type

Delphi仿速达财务软件导航条组件开发教程

Delphi作为一款历史悠久的集成开发环境(IDE),由Embarcadero Technologies公司开发,它使用Object Pascal语言,被广泛应用于Windows平台下的桌面应用程序开发。在Delphi中开发组件是一项核心技术,它允许开发者创建可复用的代码单元,提高开发效率和软件模块化水平。本文将详细介绍如何在Delphi环境下仿制速达财务软件中的导航条组件,这不仅涉及到组件的创建和使用,还会涉及界面设计和事件处理等技术点。 首先,需要了解Delphi组件的基本概念。在Delphi中,组件是一种特殊的对象,它们被放置在窗体(Form)上,可以响应用户操作并进行交互。组件可以是可视的,也可以是不可视的,可视组件在设计时就能在窗体上看到,如按钮、编辑框等;不可视组件则主要用于后台服务,如定时器、数据库连接等。组件的源码可以分为接口部分和实现部分,接口部分描述组件的属性和方法,实现部分包含方法的具体代码。 在开发仿速达财务软件的导航条组件时,我们需要关注以下几个方面的知识点: 1. 组件的继承体系 仿制组件首先需要确定继承体系。在Delphi中,大多数可视组件都继承自TControl或其子类,如TPanel、TButton等。导航条组件通常会继承自TPanel或者TWinControl,这取决于导航条是否需要支持子组件的放置。如果导航条只是单纯的一个显示区域,TPanel即可满足需求;如果导航条上有多个按钮或其他控件,可能需要继承自TWinControl以提供对子组件的支持。 2. 界面设计与绘制 组件的外观和交互是用户的第一印象。在Delphi中,可视组件的界面主要通过重写OnPaint事件来完成。Delphi提供了丰富的绘图工具,如Canvas对象,使用它可以绘制各种图形,如直线、矩形、椭圆等,并且可以对字体、颜色进行设置。对于导航条,可能需要绘制背景图案、分隔线条、选中状态的高亮等。 3. 事件处理 导航条组件需要响应用户的交互操作,例如鼠标点击事件。在Delphi中,可以通过重写组件的OnClick事件来响应用户的点击操作,进而实现导航条的导航功能。如果导航条上的项目较多,还可能需要考虑使用滚动条,让更多的导航项能够显示在窗体上。 4. 用户自定义属性和方法 为了使组件更加灵活和强大,开发者通常会为组件添加自定义的属性和方法。在导航条组件中,开发者可能会添加属性来定义按钮个数、按钮文本、按钮位置等;同时可能会添加方法来处理特定的事件,如自动调整按钮位置以适应不同的显示尺寸等。 5. 数据绑定和状态同步 在财务软件中,导航条往往需要与软件其他部分的状态进行同步。例如,用户当前所处的功能模块会影响导航条上相应项目的选中状态。这通常涉及到数据绑定技术,Delphi支持组件间的属性绑定,通过数据绑定可以轻松实现组件状态的同步。 6. 导航条组件的封装和发布 开发完毕后,组件需要被封装成独立的单元供其他项目使用。封装通常涉及将组件源码保存为pas文件,并在设计时能够在组件面板中找到。发布组件可能还需要编写相应的安装包和使用文档,方便其他开发者安装和使用。 7. Delphi IDE的支持 Delphi IDE提供了组件面板编辑器(Component Palette),允许开发者将开发好的组件添加到组件面板中。在组件面板编辑器中,可以自定义组件的图标和分类,使得组件在Delphi中的使用更为便捷。 通过以上的知识点梳理,可以看出Delphi仿速达导航条组件的开发涉及到的不仅仅是简单的代码编写,还涉及到用户界面设计、事件驱动编程、组件封装等多个方面。掌握这些知识点,对于一名Delphi开发者而言,是十分重要的。
recommend-type

【空间分布规律】:甘肃土壤类型与农业生产的关联性研究

# 摘要 本文对甘肃土壤类型及其在农业生产中的作用进行了系统性研究。首先概述了甘肃土壤类型的基础理论,并探讨了土壤类型与农业生产的理论联系。通过GIS技术分析,本文详细阐述了甘肃土壤的空间分布规律,并对其特征和影响因素进行了深入分析。此外,本文还研究了甘肃土壤类型对农业生产实际影响,包括不同区域土壤改良和作物种植案例,以及土壤养分、水分管理对作物生长周期和产量的具体影响。最后,提出了促进甘肃土壤与农业可持续发展的策略,包括土壤保护、退化防治对策以及土壤类型优化与农业创新的结合。本文旨在为