前回から続いて、読み込んだcsvファイルから薬剤内服の有無、副作用ADRの有無、を抽出して、disproportionality analysisを行うまでの一例を記載する。


１。Disproportionality analysisとは

検索すれば色々な資料で説明がかなりなされているので簡単に。

普通のコホートデータなどで、薬と症状の関係を調べる場合、基本的には、
薬内服の有無 VS 症状発現の有無　で２＊２分割表を作って、カイ2乗検定やフィッシャー正確検定をするのが一般的である。
（表１）

適当に作った上の表１だと、フィッシャー正確検定では

> fisher.test(matrix(c(30,300,300,5000),ncol=2))
	Fisher's Exact Test for Count Data

data:  matrix(c(30, 300, 300, 5000), ncol = 2)
p-value = 0.01503
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 1.085286 2.480727
sample estimates:
odds ratio 
  1.666477 

 のようになるので、Odds比 1.67 （95％信頼区間 1.085~2.481）で下95%>1なので、薬服用と症状発現には有意な関連があることになる。

ただし、JADERなどの自発報告データでは、基本的には薬を内服して副作用が出た人だけ報告されるため、表１の左上のマスに該当する症例しか報告されず、普通には解析ができない。そのため、特定の薬剤Aを服用したかどうか＊特定の症状Bが出たかどうか、の分割表にして同じように検討することになる。

これがdisproportionality analysisで、この場合はORではなくReporting OR (=ROR)と呼ばれることになる。
ORの場合は上95%<1の場合は有意に低いOR、下95%>1の場合は有意に高いOR、という解釈ができるが、RORの場合は有意に高いかどうかだけが着目されている。
ROR >1のときの解釈として、「その症状ADRはその薬の服用に伴って有意に（他の様々な薬による場合よりも）多く報告されている（= highly reported）」ということになる。つまり、そのような薬剤は、その症状の発現リスクが実際に高いかもしれないので、今後コホート研究などで確かめるべき下地にできる、ということになる。（逆にROR<1の時に、その症状はその薬剤によってリスクが下がるかもしれない？、という解釈は一般的ではないと思われる）
ここで、内服しなかった場合と比べているわけではないということ、そのため、その症状が起こるfrequencyはわからないこと（服用した全体の分母がわからないので当然）、またある薬AによるRORとある薬BによるRORとで多寡を薬剤間で比較することにあまり大きな意味づけはできないということ、には注意が必要である。



２。具体例

次に、具体的にある薬剤とある症状の関係を簡便にチェックしたい場合のRでの操作を述べる。基本的には、報告ごとに服用の有無、また症状の有無について0/1の2値ベクトルでラベルして、分割表にすれば良い。問題となるのが、数十万単位以上の件数の各報告に対して服用・症状の抽出方法で、forループでやると絶望的に遅くなるため、可能な限り、ベクトル操作ないしapply利用でやる必要がある。


例えば、スルピリドと、その内服により高頻に起きることが既知の薬剤性パーキンソニズムについて調べてみたいとする。重複は最低限除外しているとして、

> tab <- as.data.frame(table(drug[["医薬品.一般名."]]))

とすれば薬剤一般名の内訳がわかる。「スルピリド」を含むものは

 > drug_t <-  drug[ grep("スルピリド", drug[["医薬品.一般名."]]) , ] ; table(drug_t[,"医薬品.一般名."])

スルピリド 
      4857 

と1種類しかないので、 

> demo[ intersect(demo$ID,  drug[ grep( "スルピリド" , drug[["医薬品.一般名."]] ) , "ID"]),  "sulpiride"] <- 1 
> demo$sulpiride <- ifelse( is.na(demo$sulpiride )==T, 0, 1)
> table(demo$sulpiride)

     0      1 
639800   4857 

のように症例単位でスルピリド使用歴の有無（0/1）で2値化できる。 
のようになるので、「パーキンソニズム・異常高熱」とかぶるため、「パーキンソニズム」だけを検索できるようにすると良さそう、となる。そこで、

> demo[ intersect(demo$ID,  reac[ reac[["有害事象"]] %in% c("パーキンソニズム") , "ID"]) ,  "parkinsonism"] <- 1    #
> demo$parkinsonism <- ifelse( is.na(demo$parkinsonism)==T, 0, 1)
> table(demo$parkinsonism)

     0      1 
643899     758
 

とすると、パーキンソニズムを発現したかどうかを各症例にラベリングできる。
このまま分割表にすると

> table(demo$sulpiride, demo$parkinsonism)
   
         0      1
  0 639213    587
  1   4686    171

となり、さらにRORを計算すると
 

> fisher.test(table(demo$sulpiride, demo$parkinsonism))

	Fisher's Exact Test for Count Data

data:  table(demo$sulpiride, demo$parkinsonism)
p-value < 2.2e-16
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 33.23865 47.25546
sample estimates:
odds ratio 
  39.71508 

で ROR 39.7（95% CI: 33.2 ~ 47.3）  となるので、スルピリドに伴ってパーキンソニズムは有意に多く報告されている（significantly highly reported）ということになる。これは既報告と整合する結果と言える。
なお上記のコードだと投与日→投与日という経過が必ずしも保証されていないことには注意が必要である。

「パーキンソニズム」を検索するところのコードで %in% を使ったがこれにはdplyrが必要。 %in% c( "A", "B")と置けば　AないしBの副作用に当てはまるものをラベルすることになるので、副作用カテゴリーの有無、というラベルの仕方もできる（パーキンソニズム関連：「パーキンソニズム」ないし「安静時振戦」、など）。

また薬剤を指定するところのコードではgrep（"パーキンソ", ...）を使っているが、grep（"A|B", ...）などとすればAないしBの文字列に一致するものを検索できる。前述のように特異性には気をつける必要がある。病名と同様に、%in%で直接指定する方が予期せぬ間違いをしてしまうという恐れは少ない。