我现在收集了两组不同类型的数据,想像在cma软件中一样,输入这两组不同的数据,然后合起来进行meta分析,我采用了下面的代码
# 加载必要的包
library(readxl)
library(metafor)
library(dplyr) # 用于简化数据处理
library(knitr)
library(htmltools)
library(kableExtra)
library(forestploter)
### 1. 处理 Mean 数据
# 读取 mean 文件
data_means <- read_excel(file.choose()) # 选择 mean 文件
# 提取均值、标准差和样本量
m1 <- data_means$`Group-A mean`
sd1 <- data_means$`Group-A Std-Dev`
n1 <- data_means$`Group-A sample size`
m2 <- data_means$`Group-B mean`
sd2 <- data_means$`Group-B Std-Dev`
n2 <- data_means$`Group-B sample size`
# 计算自由度和校正因子
df_means <- n1 + n2 - 2
J_means <- 1 - (3 / (4 * df_means - 1))
# 计算 Cohen's d 和校正后的 Hedges' g
SMD <- (m1 - m2) / sqrt(((n1 - 1) * sd1^2 + (n2 - 1) * sd2^2) / df_means)
hedges_g_means <- SMD * J_means
# 计算标准误和方差
stderr_d_means <- sqrt((n1 + n2) / (n1 * n2) + (SMD^2 / (2 * df_means)))
stderr_g_means <- stderr_d_means * J_means
variance_g_means <- stderr_g_means^2
# 创建 mean 数据框
dat_means <- data.frame(
yi = hedges_g_means,
StdErr = stderr_g_means,
Variance = variance_g_means,
stringsAsFactors = FALSE
) %>%
cbind(data_means %>% select(-starts_with("Group"))) # 保留其他列
### 2. 处理 Correlation 数据
# 读取 correlation 文件
data_corr <- read_excel(file.choose()) # 选择 correlation 文件
# 提取相关系数和样本量
r <- data_corr$correlation
N <- data_corr$`sample size`
# 计算 Fisher's Z 变换效应量
dat_corr <- escalc(measure="ZCOR", ri=r, ni=N)
# 创建 correlation 数据框
dat_corr <- data.frame(
yi = dat_corr$yi,
StdErr = sqrt(dat_corr$vi), # 标准误
Variance = dat_corr$vi, # 方差
stringsAsFactors = FALSE
) %>%
cbind(data_corr %>% select(-correlation, -`sample size`))
### 3. 合并两组数据
# 找到两个数据框的公共列
common_cols <- intersect(names(dat_means), names(dat_corr))
# 筛选公共列,确保结构一致
dat_means_filtered <- dat_means %>% select(all_of(common_cols))
dat_corr_filtered <- dat_corr %>% select(all_of(common_cols))
# 合并两组数据
dat_combined <- rbind(dat_means_filtered, dat_corr_filtered)
# 移动 'study name' 到第一列
dat_combined <- dat_combined %>%
select(`Study name`, everything()) # 将 'study name' 列移到第一列
# 进行随机效应模型的Meta分析
res_combined <- rma(yi, Variance, data = dat_combined, method = 'REML')
最终得到的随机效应模型的效应量的值与cma相比相差很大,如果单独只使用一个类型的数据测试,则不会出现大的偏差,想知道我该怎么修改r的代码,才能正确进行meta分析