原文連結：http://tecdat.cn/?p=24861

原文出處：拓端資料部落公眾號

概括

1. r 語言允許使用者計算 lme 4 包中廣義線性混合模型的功效。功率計算基於蒙特卡羅模擬。
2. 它包括用於 (i) 對給定模型和設計進行功效分析的工具；(ii) 計算功效曲線以評估功效和樣本量之間的權衡。
3. 本文提供了一個教程，使用具有混合效果的計數資料的簡單示例（具有代表環境監測資料的結構）。

介紹

假設檢驗的功效定義為假設原假設為假，檢驗拒絕原假設的概率。換句話說，如果一個效應是真實的，那麼分析判斷該效應具有統計顯著性的概率是多少？

如果一項研究的功效不足，資源可能被浪費，真正的效果可能被遺漏。另一方面，一項大型研究的花費可能過大，因此其費用也會超過必要的範圍。因此，在收集資料之前進行功效分析是一個很好的做法，以確保樣本具有適當的規模來回答正在考慮的任何研究問題。

廣義線性混合模型 (GLMM) 在生態學中很重要，它允許分析計數和比例以及連續資料，並控制空間非獨立性.

蒙特卡羅模擬是一種靈活且準確的方法，適用於現實的生態研究設計。在某些情況下，我們可以使用解析公式來計算功效，但這些通常是近似值或需要特殊形式的設計。模擬是一種適用於各種模型和方法的單一方法。即使公式可用於特定模型和設計，定位和應用適當的公式也可能非常困難，因此首選模擬。

對於對 r 不夠熟悉的研究人員，設定模擬實驗可能太複雜了。在本文中，我們介紹了一個工具來自動化這個過程。

r 包

有一系列的 r 包目前可用於混合模型的功效分析 。然而，沒有一個可以同時處理非正態因變數和廣泛的固定和隨機效應規範。

圖1

r 旨在與任何可以與 lme 4 中的 lmer 或 glmer 配合的線性混合模型 (LMM) 或 GLMM 一起使用。這允許具有不同固定和隨機效應規範的各種模型。還支援在 r 中使用 lm 和 glm 的線性模型和廣義線性模型，以允許沒有隨機效應的模型。

r 中的功效分析從適合 lme 4 的模型開始。

在 r 中，通過重複以下三個步驟來計算功效：(i) 使用提供的模型模擬因變數的新值；(ii) 將模型重新擬合為模擬因變數；(iii) 對模擬擬合應用統計檢驗。在此設定中，已知存在測試效果，因此每個陽性測試都是真正的陽性，每個陰性測試都是 II 類錯誤。可以根據步驟 3 的成功和失敗次數計算測試的功效。

教程

本教程使用包含的資料集。該資料集代表環境監測資料，在連續固定效應變數x（例如研究年份）的10 個水平上測量三個組g（例如研究地點）的因變數z（例如鳥類丰度）。還有一個連續因變數y，在本教程中沒有使用。

擬合模型

我們首先將 lme 4 中的一個非常簡單的泊松混合效應模型擬合到資料集。在這種情況下，我們有一個隨機截距模型，其中每個組 (g) 都有自己的截距，但這些組共享一個共同的趨勢。

1. glm
2. summary

本教程重點介紹關於x趨勢的推斷。在這種情況下，x的估計效應大小為-0.11，使用預設z檢驗在 0.01 水平上顯著。

請注意，我們特意使用了一個非常簡單的模型來使本文易於理解。例如，適當的分析會包含更多的組，並會考慮過度分散等問題。。

簡單的功率分析

假設我們想重複這項研究。如果效果是真實的，我們是否有足夠的功效來期待積極的結果？

指定效應量

在開始功效分析之前，重要的是要考慮您感興趣的效果大小型別。功效通常隨效果大小而增加，較大的效果更容易檢測。回顧性“觀察功效”計算，其中目標效應大小來自資料，給出誤導性結果.

對於此示例，我們將考慮檢測 -0.05 斜率的功效。可以使用 lme 4 函式擬合 glmer 模型中的固定效應。然後可以更改固定效應的大小。變數x的固定效應的大小可以從 -0.11 更改為 -0.05，如下所示：

3. fixe<‐ ‐0.05

在本教程中，我們只更改變數x的固定斜率。但是，我們也可以更改隨機效應引數或殘差方差（適用於合適的模型）。

執行功效分析

一旦指定了模型和效應大小，在 r 中進行功效分析就非常容易了。由於這些計算基於蒙特卡羅模擬，因此您的結果可能略有不同。如果你想得到和教程一樣的結果，你可以使用 set.seed(123)。

1. power

鑑於此特定設定，拒絕x中零趨勢的零假設的能力約為 33%。這幾乎總是被認為是不夠的；傳統上，80% 的功率被認為是足夠的.

在實踐中，z檢驗可能不適合這樣一個小例子。引數引導測試可能是最終分析的首選。但是，更快的z-test 更適合學習使用該包以及在功效分析期間進行初始探索性​​工作。

增加樣本量

在第一個示例中，估計功率很低。小型試點研究通常沒有足夠的功效來檢測微小的影響，但更大的研究可能會。

試點研究對x 的10 個值進行了觀察，例如代表研究第 1 年到第 10 年。在此步驟中，我們將計算將其增加到 20 年的影響。

1. modl2 <‐ extend
2. power(modl2)

沿引數指定要擴充套件的變數，n 指定要替換它的級別。擴充套件模型 2 現在將具有從 1 到 20 的x值，與以前一樣分為三組，總共 60 行（與模型 1 中的 30 行相比）。

通過觀察x 的20 個值，我們將有足夠的能力來檢測大小為 -0.05 的效應。

各種樣本量的功效分析

當資料收整合本高昂時，使用者可能只想收集達到一定統計能力所需的資料量。 功效曲線 函式可用於探索樣本大小和功效之間的權衡。

確定所需的最小樣本量

在前面的示例中，當對變數x 的20 個值進行觀察時，我們發現了非常高的功效。我們能否減少這個數字，同時保持我們的功效高於通常的 80% 閾值？

1. poerCure
3. print
4. plot

請注意，我們已將此結果儲存到變數 pc2 以匹配模型 2 中的編號。由於模型 1 沒有足夠的功率，我們沒有通過 powerCurve 執行它。繪製的輸出如圖所示。我們可以看到，檢測x趨勢的能力隨著取樣大小的增加而增加。這裡的結果基於將模型擬合到 10 個不同的自動選擇的子集。最小的子集僅使用前 3 年（即 9 個觀測值），最大的子集使用所有 20 個假設研究年份（即 60 行資料）。該分析表明，該研究必須執行 16 年才能有≥80% 的功效來檢測指定大小的影響。

圖2

檢測大小為 -0.05 的固定效應的功效 (±95% CI)，使用 powerCurve 函式在一系列樣本大小上計算。變數x的不同值的數量從 3 (n= 9) 到 20 (n= 60) 不等。

改變組的數量和大小

增加觀察到的x值的數量可能不可行。例如，如果x是研究年份，我們可能不願意等待更長時間的結果。在這種情況下，增加研究地點的數量或每個地點的測量數量可能是更好的選擇。這兩項分析從我們的原始模型 1 開始，該模型已有 10 年的研究時間。

新增更多組

我們可以像為x新增額外值一樣為g新增額外級別。例如，如果變數g代表我們的研究站點，我們可以將站點數量從 3 增加到 15。

1. extend(n=15)
2. plot(pc3)

與上一個示例的主要變化是我們將變數g傳遞給了沿引數。該分析的輸出如圖 1 所示。要達到 80% 的功率，我們至少需要 11 個站點。

圖 3

檢測大小為 -0.05 的固定效應的功效 (±95% CI)，使用 powerCurve 在一系列樣本大小上計算。因子g的級別數從 3 (n= 30) 到 15 (n= 150) 不等。

增加組內的大小

我們可以用內參數替換擴充套件和 powerCurve 的沿引數以增加組內的樣本大小。每個組在x和g 的每個水平上只有一個觀察值。我們可以將其擴充套件到每個站點每年 5 次觀測，​​如下所示：

1. extend( n=5)
3. plot(p4)

請注意 powerCurve 的breaks 引數。為x和g 的每個組合提供一到五個觀察結果。圖表明每年每個站點 4 次觀測會給我們 80% 的效力。

圖 4

檢測大小為 -0.05 的固定效應的功效 (±95% CI)，使用 powerCurve 函式在一系列樣本大小上計算。x和g 的每個組合的觀察數從 1 (n= 30) 到 5 (n= 150) 不等。

最受歡迎的見解

1.Matlab馬爾可夫鏈蒙特卡羅法（MCMC）估計隨機波動率（SV，Stochastic Volatility） 模型

2.基於R語言的疾病製圖中自適應核密度估計的閾值選擇方法

3.WinBUGS對多元隨機波動率模型：貝葉斯估計與模型比較

4.R語言迴歸中的hosmer-lemeshow擬合優度檢驗

5.matlab實現MCMC的馬爾可夫切換ARMA – GARCH模型估計

6.R語言區間資料迴歸分析

7.R語言WALD檢驗 VS 似然比檢驗

8.python用線性迴歸預測股票價格

9.R語言如何在生存分析與Cox迴歸中計算IDI，NRI指標