發表文章

目前顯示的是 5月, 2024的文章

假設檢定

虛無假設:效果量=0 對立假設:效果量≠0 p 值

貝氏分析

貝氏定理:P(A|B) = P(A,B) / P(B) = P(B|A) * P(A)/P(B), p(θ|D)= p(θ, D)/p(D) = p(θ) x p(D|θ)/p(D), theta: 母數, D:資料或變數 x, 條件機率 = 聯合機率/邊際機率 邊際機率:P(A), P(B), p(D), 不管/積分掉其他的變數 條件機率:P(A|B), P(B|A), p(θ|D), p(D|θ) 聯合機率:P(A,B) = P(A|B) * P(B) = P(B|A) * P(A) 如果兩個事件 A 和 B 是獨立的,那麼條件機率等於邊際機率 似然函數 L(θ) ∝ ∏p(xi|θ), L(θ|x) ∝ p(x|θ):看見 x 時,θ 的似然率等比於 θ 的分布(例如平均 0、標準差 1 的標準常態 Z 分布)固定時 x 的密度、 機率 vs 似然率:單一宇宙 vs 多維宇宙 機率密度函數 f(x) or pdf:dnorm(x, mean, sd), p(x|θ), 積分等於 1 累積分布函數 cdf:pnorm(x, mean, sd) 分位數函數(pnorm 的 反函數):qnorm(0.975) 頻率學派 vs 貝氏學派:機率是重複事件發生的頻率 vs 主觀上某事件發生的不確定性、似然率 L(θ|x) 等比於 p(D|θ) vs p(D|θ)、樣本是隨機的變數 vs 已經發生的固定數、母數(參數)是未知的常數(沒有機率可言) vs 未知的隨機變數(有機率)、θ 固定 vs p(θ)、客觀(只依據本次試驗的結果)vs 主觀(依據先驗信念及本次試驗的結果而更新成後驗信念)、傳統 vs 新、簡單vs 複雜、人工計算可行 vs 不可行、最大似然率估計 vs 邊際似然率 p(D) = ∫p(θ’)p(D|θ’)dθ’ 離散機率的總和等於 1、連續機率的積分等於 1、似然率的總和/積分不等於 1。 邊際似然率是所有可能 θ 的積分而非所有似然率的平均值。 H0: θ=0, H1: θ ≠ 0 In fact, p(θ|D) is p(θ|D, Mi), p(θ) is p(θ|Mi), p(D|θ) is p(D|θ, Mi), p(D) is p(D|Mi) Under H1 or f1(x|θ1): p(H1|x) = p(H1) x p(x|H1)/p(x)...

頻率學派 vs 貝氏學派

機率是重複事件長期發生的頻率 vs 主觀上某事件發生的不確定性(比較像是一般人的觀念)。 傳統 vs 新 簡單vs 複雜 人工計算可行 vs 不可行 母數(參數)是未知的常數(沒有機率可言) vs 未知的隨機變數(有機率) 樣本是隨機的變數 vs 已經發生的固定數 客觀(只依據本次試驗的結果)vs 主觀(依據先驗信念及本次試驗的結果而更新成後驗信念)。 H: 假設,D: 資料。H0: 虛無假設(參數或母數 θ 等於 0),H1: 對立假設,p(H1)=1-p(H0)。 一個棒球比賽的觀眾關心的是「機率」p(D|H),亦即預測「這個左打者擊中球的機率是多少?」;當他看見一名打者擊中了球以後,計算這一名打者是一個左打者的機率是多少則是屬於「逆機率」p(H|D)。 醫生關心的是「逆機率」p(H|D):當看見某症候、檢驗、檢查 D時,用貝氏定理計算病人有某疾病 H 的機率。 估計 p 值(H0 為真時,看見 D 或比它更極端的機率 p(D|H0)) vs 估計後驗機率 p(H1|D):p(D|H) ≠ p(H|D)。 95% 信賴區間包含母數的機率是 0% 或 100% vs 95% 置信區間包含母數的機率是 95%。 頻率學派:當 p < α(例如 0.05)時拒絕 H0,當 p ≥ α 時無法拒絕 H0。但是我們無法知道 p(H0|D)、p(H1|D),而且無論 p 值多大都無法接受 H0。 貝氏定理: 1. H1 的後驗勝算 p(H1|D)/(p(H0|D))= 先驗勝算 x 似然比或貝氏因子 p(D|H1)/p(D|H0),勝算=機率/(1-機率),機率=勝算/(勝算+1)。 2. 後驗機率 p(H|D)= 先驗機率 p(H) x p(D/H)/p(D) p(D): p(D/H0) + p(D/H1) RCT 時,頻率學派的 α 是 H0 為真(無效)時,你卻斷言為假(有效)的機率:p(斷言有效|無效)。傳統上它被稱為「第一型錯誤率」,其實它應該是「第一型斷言率」,而不是治療無效你卻以為有效的錯誤機率:貝氏學派的 1 - p(有效|資料,先驗機率)。後者才是我們所關心的。 ESP and the significance of significance 22 https://understandinguncertainty.org/node/1286 爲什麼我...

ACCEPT analysis

# ACCEPT: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7613267/ source(url("https://github.com/guhjy/guhjy/raw/main/shrinkage.R")) # ORBITA: https://www.semanticscholar.org/paper/Percutaneous-coronary-intervention-in-stable-angina-Al-Lamee-Thompson/1ba069c9f4a819503735daa7ca8abdcbd69f9200 onetailed_p <- function(estimate, lower, upper) {   se <- (upper - lower) / 3.92   z <- estimate/se;p <- pnorm(z)   cat("One-tailed p-value (more than):", p, "\n")   return(p)} onetailed_p(16.6, -8.9, 42) # DanGer Shock: https://www.nejm.org/doi/full/10.1056/NEJMoa2312572 onetailed_p <- function(estimate, lower, upper) {    se <- (log(upper) - log(lower)) / 3.92    z <- log(estimate)/se;p <- pnorm(z)    cat("One-tailed p-value (more than):", p, "\n")    return(p)} onetailed_p(0.74, 0.55, 0.99) accept <- function(estimate, lower_ci, upper_ci, threshold) {   se <- (log(upper_ci) - log(lower_ci)) / (2 * qnorm(0....

Shrinkage in Claude 3 Opus

# 連續變項時,在 Mistral-large 上載此 pdf: Evaluating a shrinkage estimator for the treatment effect in clinical trials 23 https://onlinelibrary.wiley.com/doi/10.1002/sim.9992 用中文只顯示收縮估計及其標準誤與信賴區間: b=0.4, s=0.3 # Percutaneous coronary intervention in stable angina (ORBITA): a double-blind, randomised controlled trial 17 用中文只顯示收縮估計及其標準誤與信賴區間: b=16.6, s=13 GPT 4: RCT 的 estimate 16.6 (95% confidence interval -8.9, 42),使用 pymc3 進行貝氏再分析,在弱信息先驗條件下抽樣 500000 次,計算後驗分布及 estimate 小於 0 的概率,並進行敏感性分析,只用繁體中文顯示結果,而不顯示代碼 # HR, OR, RR 時,在 Claude 3 Opus 上載此 pdf:  Microaxial Flow Pump or Standard Care in Infarct-Related Cardiogenic Shock 24 https://www.nejm.org/doi/full/10.1056/NEJMoa2312572 Convert hr 0.74 (95% ci 0.55, 0.99) to ln(hr) as the b and se of ln(hr) as the s to calculate shrinkage est and 95% ci of hr and show the results without showing the procedures and recheck the results RCT 的 hr 0.74 (95% confidence interval 0.55, 0.99),使用 pymc3 進行貝氏再分析,在弱信息先驗條件下抽樣 500000 次,計算後驗分布及 hr 小於 1 的概率,並進行敏感...

Shrinkage estimator datacamp

source(url("https://github.com/guhjy/guhjy/raw/main/shrinkage.R")) result <- shrinkage_estimator(b = 0.4, s = 0.3) cat("收縮估計:", result$betahat, "\n") cat("收縮估計的標準誤:", result$sigma_hat, "\n") cat("95% 信賴區間: [", result$confidence_interval[1], ", ", result$confidence_interval[2], "]\n")

Shrinkage estimator

 shrinkage_estimator <- function(b, s) {   z <- b / s   p <- c(0.32, 0.31, 0.3, 0.07)   sigma <- c(0.61, 1.42, 2.16, 5.64)   sigma2 <- sigma^2   # mixing probabilities   q <- p * dnorm(z, 0, sqrt(sigma2 + 1))   q <- q / sum(q)  # Normalize to sum to 1   # conditional means and variances   pm <- b * sigma2 / (sigma2 + 1)   pv <- s^2 * sigma2 / (sigma2 + 1)   # shrinkage estimator (β̂)   betahat <- sum(q * pm)   # second moment and total variance for the mixture   pm2 <- sum(q * pm^2)   ps2 <- sum(q * pv)   # standard deviation (σ) of the mixture distribution   sigma_hat <- sqrt(ps2 + pm2 - betahat^2)   # 95% confidence interval around β̂   ci_lower <- betahat - 1.96 * sigma_hat   ci_upper <- betahat + 1.96 * sigma_hat   # Return results as a list   list(     betahat = betahat,     sigma_hat = sigma_hat, ...

Specification curve analysis

#  Grilling the data: Application of specification curve analysis to red meat and all-cause mortality 24 : r #  specr::setup install.packages("specr") install.packages("tidyverse") library(tidyverse) library(specr) # Function to create logistic regression models log_glm <- function(formula, data) {   glm(formula, data = data, family = binomial()) } # Simulate example data set.seed(123)  # For reproducibility n <- 1000  # Number of observations example_data <- tibble(   x1 = rnorm(n, mean = 50, sd = 10),  # Simulate red meat consumption   c1 = sample(c("Male", "Female"), n, replace = TRUE),  # Sex   c2 = rnorm(n, mean = 50, sd = 10),  # Age   c3 = rnorm(n, mean = 25, sd = 5),   # BMI   c4 = rbinom(n, 1, 0.2),             # Diabetes status (0 = no, 1 = yes)   y1 = as.factor(rbinom(n, 1, plogis(-3 + 0.02 * x1 + 0.1 * c4)))  # Mortality (logistic model) ) # Setup the specif...