
---
title: "整合分析（Meta-Analysis）"
author: "許清芳，鄭中平"
date: 2024-04-27
output: slidy_presentation
---

## 小題

(1)請用下列方式讀入資料，並繪製森林圖。
dta <- read.csv("https://raw.githubusercontent.com/ChungPingCheng/R4BS2/main/Data/Deindividuation.csv", header=T) 

(2)請計算平均效果量，並檢驗是否為0。

(3)請繪製漏斗圖，並檢驗是否有出版偏誤。

(4)請檢驗匿名團體的參與者是否在群體中，是否會影響研究中效果的大小。

### 提示

先找到範例中對應的程式複製到新的 RMD，再往下修改

## 看一下資料
```{r pkgs, echo=TRUE, cache=FALSE, message=FALSE, warning=FALSE}
#整體設定，含載入套件
source("https://raw.githubusercontent.com/ChungPingCheng/R4BS2/main/R4BS_setup.R")
```


```{r }
#讀進資料
dta <- read.csv("https://raw.githubusercontent.com/ChungPingCheng/R4BS2/main/Data/Deindividuation.csv", header=T) 
```


```{r }
#顯示前六筆資料
dim(dta)
head(dta)
```

## (1)請用下列方式讀入資料，並繪製森林圖。

### 描述性統計

```{r eval=FALSE}
#依不同Par_in_groups，看看效果量等的描述性統計
dta |>
  select(-Study) |>
  gtsummary::tbl_summary(by=Par_in_groups,
                statistic = list(all_continuous() ~ "{mean} ({sd})"))
```




### 資料圖形化

```{r}
#效果量直方圖
ggplot(data=dta, 
       aes(r, after_stat(density)))+
  geom_histogram(binwidth = function(x) 2 * IQR(x) / (length(x)^(1/3)),
                 col=1, fill=8)+
  labs(x="研究所得效果（相關係數量）",
       y="密度")
```

## (2)請計算平均效果量，並檢驗是否為0。

### 整合分析 

```{r }
#利用相關與樣本數計算相關的抽樣變異，並放進原先資料
dta <- metafor::escalc(measure="COR", ri=r, ni=N, data=dta)
head(dta)
```

```{r }
#計算所有研究的平均效果
summary(res <- metafor::rma(yi=yi, vi=vi, data=dta))
```

### 第一小題
### 森林圖，圖示結果
```{r}
metaviz::viz_forest(res)
```

## (3)請繪製漏斗圖，並檢驗是否有出版偏誤。

### 檢查出版偏誤

### 漏斗圖 
```{r }
#檢查出版偏誤，先繪製 funnelplot
metaviz::viz_funnel(res)
```


```{r }
ranktest(res)
```

```{r }
regtest(res)
```

## (4)請檢驗匿名團體的參與者是否在群體中，是否會影響研究中效果的大小。

### 整合調節分析

```{r }
dta <- dta |>
  dplyr::mutate(Par_in_groups = factor(Par_in_groups))
```

```{r }
set.seed(201408)
dta_13 <- dta |> 
  dplyr::sample_frac(.33) |> 
  arrange(Par_in_groups, r)
```


```{r }
#圖11.5
with(dta_13, metafor::forest(yi, vi, slab = Par_in_groups))
```


```{r }
dta <- 
  dta |> arrange(Par_in_groups, r)
res_rgn <- rma(yi=yi, vi=vi, mods = ~ Par_in_groups, data=dta)
res_rgn |> summary()
```

```{r }
#繪製 forest
p1 <- metaviz::viz_forest(res_rgn, group=dta[,'Par_in_groups'],
              xlab = "r", x_limit=c(-.80,1.1), type="study_only")
p2 <- metaviz::viz_forest(res_rgn, group=dta[,'Par_in_groups'],
              xlab = "r", x_limit=c(-.80,1.1),type="summary_only")
p1/p2+plot_layout(heights = c(7, 1))
```





