ACP Japan Chapter のRFCの活動報告をシェアします
Case Crossover Design – (3)
Case Crossover Design – (3)
前回の記事のつづきです
Logistic regression model
このアプローチでは、観察期間を単位期間ごと(この場合は1hr)に分割して、イベントが発生した時以外は、コントロールとして、多くのcontrolと少数回のcaseを観察したcase-control studyのようにみなします。原文では[the case-crossover design can be viewed as matched case-control design with 1:M matched pairs]となっていまして、この、「みなす」というのが感覚的によくわからないのですがそういうことらしいです。例題では一人につき1年間観察して、1回AMIが発症たことにしていますので、1例の症例が、集計上8759例のcontrolと1例のcaseのように扱われます。例題では10例いますのでトータルの集計上の観察症例数は87600例と表示されます。matched pairを示す指標としてcase idを使用した、conditional logistic regression modelを用います。Formulaには[case~exposure+strata(id)]を入れます。
ちなみに、このスクリプトで使用されていますreshape2というライブラリのmeltは、こんなデータの変形をしたいと思うときに手作業でやってるようなデータの変形をやってくれる興味深い関数です。
matrix<-matrix(round(c(rr,lo,hi,or,lo.or,hi.or),2),nrow=2,byrow=TRUE)
rownames(matrix)<-c(“RR”,”OR”)
colnames(matrix)<-c(“Value”,”low 95% CI”,”high 95% CI”)
matrix
mat<-matrix(, nrow = T-1, ncol = 0)
for (i in 1:10) {
if (T%%frq[i]==0) {
exposure<-c(rep(c(1,rep(0,T/frq[i]-1)),frq[i]))[-T]
} else {
exposure<-c(rep(c(1,rep(0,trunc(T/frq[i])-
1)),frq[i]),rep(0,T-frq[i]*trunc(T/frq[i])))[-T]
}
mat<-cbind(mat,exposure)
}
commat<-rbind(efftime,mat)
library(reshape2)
data.wide<-as.data.frame(commat)
colnames(data.wide)<-c(1:10)
data<-melt(data.wide,measure=c(1:10))
colnames(data)<-c(“id”,”exposure”)
data$case<-rep(c(1,rep(0,T-1)),5)
head(data)
library(survival)
mod<- clogit(case~exposure+strata(id),data)
summary(mod)
結果として得られたodds ratioは、前記事の未調整のORと大体同じになります。
Time trend adjustment with conditional logistic regression model
Case crossover designではcase自身をcontrolとして用いるため、年齢・性別・既往歴や基礎疾患といった症例の背景因子のような交絡因子は調整されているとして扱われています。しかしながら、トレンドの交絡を避けることができません。トレンドの交絡とは?ということですが、参考にしている文献(中国論文1とMaclure文献 2)では、MIの発現は日内変動のパターンを取ることが知られていることを例にしています。そしてコーヒーを飲むのもおそらく日内変動のパターンを取りそうです。そこで、朝、昼、午後、夜に対応して1, 2, 3, 4の値を取る変数clockを導入して調整します。(ただし、例題ではAMIの発症時刻が不明ですので、データがないということでwarningが出ます。私の感覚ではvariable にmissing dataがあるとその観察は計算から除外され、自由度が減るという方がしっくりきます。)
data$clock<-c(rep(1,T/(3654)),rep(2,T/(3654)),rep(3,T/(3654)),rep(4,T/(3654)))
mod.adj<- clogit(case~exposure+clock+strata(id),data)
summary(mod.adj)
内科学会地方会2018年度スケジュール
とりあえずメモ
<http://www.naika.or.jp/jsim_wp/wp-content/uploads/2015/05/03_2018keikaku.pdf>
http://www.naika.or.jp/jsim_wp/wp-content/uploads/2015/05/03_2018keikaku.pdf
vigiRank beyond the BCPNN signal detection
先日の記事でBCPNNは、WHOが使用しているシグナル検出方法だ、というようなことを書きましたが実はWHOは少し軸足を別のものに移してきています。その名前はvigiRank。1,2
vigiRankとはどういうものでしょうか。まずは、説明です。3
vigiRank is a data-driven predictive model for emerging safety signals. In addition to disproportionate reporting patterns, it also accounts for the completeness, recency, and geographic spread of individual case reporting, as well as the availability of case narratives. Previous retrospective analysis suggested that vigiRank performed better than disproportionality analysis alone.
曰く、新興安全性シグナルを検出する、データ駆動予測モデル。不均衡な報告パターンに加えて、個別症例報告の以下の点を考慮する。
- 完全性 4(具体的には性別・年齢・被疑薬の使用理由・投与開始日・イベント発現日・narrative・dechallenge rechallenge情報や転帰情報の有無を織り込んだ指標)
- 最新性 (文献では直近3年の報告の件数(割合?)が示されていました)
- 地理的な広がり (文献では複数国からの報告が例示されていました)
- 病歴が入手されているか(定型文や断片的な情報は未入手扱い)
単純なdisproportionality分析より優れたパフォーマンスを示す可能性が示唆されている。
使うパラメータが増えた分利用できる症例数は減りそうです。次の動画の最後の方で紹介されています。この演者は全体にわかりやすいはっきりとした発音なのですが、WHOを「ダブリュエイチオー」ではなく、「フー」と発音しているので聞き落しなく。
<http://www.imi-protect.eu/documents/PROTECTSymposium-Tutorial-Statisticalsignaldetection-UMC.pdf>
References:
JT art hall
お昼休みうどんをツルんとかきこんで、JTアートホールに。
JTのビルに近づくだけでタバコくさい。中に入ったらさすがに禁煙だろうと思ったら、ビルの中もタバコを吸っている人が壁沿いにいっぱい。タバコを吸わない人間としては、死ぬかと思った。
このホールのコンサートは行かない方が良い。
BCPNN – Potential Safety Signal of Medical Product
BCPNNについて
BCPNN (A Bayesian Confidence Neural Network) 1の集計手法について概要を記します。その例として、PMDA (Pharmaceuticals and Medical Devices Agency)が公開するJADER (Japanese Adverse Drug Event Report database) データベースで、WHOが使用するシグナル検出の手法であるBCPNN (A Bayesian Confidence Neural Network) 1を計算しました。結果はこちらのページに示します。<https://wp.me/P9b6zl-gt>
BCPNNをSQLでデータ集計
SELECT Drug201711.識別番号, Drug201711.医薬品(一般名), Drug201711.医薬品の関与, Reac201711.有害事象 INTO [0010 AnalysisLineList]
FROM Drug201711 INNER JOIN Reac201711 ON Drug201711.識別番号 = Reac201711.識別番号
GROUP BY Drug201711.識別番号, Drug201711.医薬品(一般名), Drug201711.医薬品の関与, Reac201711.有害事象
HAVING (((Drug201711.医薬品の関与)=”被疑薬”));
SELECT [0010 AnalysisLineList].医薬品(一般名) AS drug, [0010 AnalysisLineList].有害事象 AS event, Count([0010 AnalysisLineList].識別番号) AS 識別番号のカウント INTO [0020 PhViDdata]
FROM [0010 AnalysisLineList]
GROUP BY [0010 AnalysisLineList].医薬品(一般名), [0010 AnalysisLineList].有害事象;
元データに入っている『”』の様な余計なごみを削除して、つぎのRスクリプトへ
Rを用いたBCPNNによるデータマイニング
setwd(“C:/Users/Oshima/Documents/2017/PMDA”)
library(PhViD)JADERdata <- read.csv(“0020PhViDdata.csv”, head=T)
PhViDdata <- as.PhViD(JADERdata)
res <- BCPNN(PhViDdata, MIN.n11=3, DECISION = 3,
RANKSTAT = 2)write.table(res$ALLSIGNALS, “BCPNN_out.csv”, quote=F, col.names=T, append=F, sep=”,”)
本質的にはこれでOK。BCPNNではない、同様の手法であるGPSでもRORも出力できますが、その辺の設定を変える場合はPhViDの取説参照。また、WHOはシグナル検出に用いる手法の軸足をBCPNNから、vigiRankという手法に移してきている模様です。
以前のs社にいたときに、これと同じ動作をするプログラムを、F社が「本当なら億の価値があるヤツだが、一緒に開発してもらうことを前提に2000万円ほどで開発しますよ」と持ちかけてきました。以前よりこんな感じのスクリプトを自分で作ってデータを見ていましたから、当然お断りしました。
Reference
Morning walk
今朝通勤中に職場近くの神社に立ち寄りました。初詣でられの準備が進んでいます。
HSP
「ヒートショック」という言葉も、ヒートショック・プロテインという物質群が報告・研究されていた時とは、全く異なる意味合いで普及してしまったな
皇居のお堀
今朝の皇居のお堀の様子です