1

私は月経周期のデータを扱っており、感染症を持っていることが月経前症状の発生を予測するかどうかを調査したいと考えています. さらに、月経前の段階が3、4、または5日(またはそれ以下またはそれ以上)続くと考える先験的な理由はありません。したがって、月経前の症状を感染の手がかりとして予測したい場合、月経前の段階に最も関連する長さを調査するために、変数「フェーズ」のさまざまなバージョンを持つモデルを比較したいと思います。

ただし、私が持っているデータと質問の種類については、AIC を計算しない glmmPQL を使用する必要があるため、MuMin や他の同様のパッケージを使用できません (qAIC も取得できません)。今のところ、ROC パッケージとパフォーマンス関数を使用してモデルを比較しましたが、それが妥当な方法であるかどうかはわかりません。以下に、データセットとモデルの詳細を示します。glmmPQLモデルを二項応答および時間的自己相関構造と比較する方法を見つけるためにネットで多くの時間を費やしましたが、私の場合は何も機能しません(たとえば、応答が二項であるためglsモデルを比較し、データが原因でlmerモデルを比較します)は自己相関しています)。どんな助けでも大歓迎ですありがとう!

データセット: 各行は 1 日を表しており、女性の中で日が繰り返されます。この例では、データは 1 月経周期のみを通知します。応答変数はバイナリ (0,1) であり、固定変数 "inf" (感染: はい/いいえ) および "段階" (段階: 月経前/その他) も同様です。データは時間的に自己相関しているため (毎日の症状は、女性の前日に発生した症状と相関しています)、時間的自己相関構造とランダム効果 ID を含めるために glmmPQL を使用します。これにより、以下のモデルが得られます (長さ: サイクルの長さ; dcycle: サイクルの日)

# model # 
library(MASS)
library(nlme)
modc.PQL<-glmmPQL(cramps~inf*phase+log(age)+log(clength), random=~1|id,  correlation=corCAR1(form=~dcycle|id), family=binomial, data=data2)

# inspecting model fit #
require(ROCR)
pr <- predict(modc.PQL)
pred <- prediction(as.vector(pr), data2$cramps)
auc.perf = performance(pred, measure = "auc")
auc.perf@y.values

ここで、さまざまなバージョンの変数 Phase を使用してモデルを比較したいと思います (phase4-> 月経の 4 日前、phase5-> 月経の 5 日前など...)。

# candidate set #
modc.PQL4<-glmmPQL(cramps~inf*phase4+log(age)+log(clength), random=~1|id, correlation=corCAR1(form=~dcycle|id), family=binomial, data=data2)      
modc.PQL5<-glmmPQL(cramps~inf*phase5+log(age)+log(clength), random=~1|id, correlation=corCAR1(form=~dcycle|id), family=binomial, data=data2)    
modc.PQL6<-glmmPQL(cramps~inf*phase6+log(age)+log(clength), random=~1|id, correlation=corCAR1(form=~dcycle|id), family=binomial, data=data2)  
modc.PQL7<-glmmPQL(cramps~inf*phase7+log(age)+log(clength), random=~1|id, correlation=corCAR1(form=~dcycle|id), family=binomial, data=data2)

# inspecting model fit model by model #
require(ROCR)
pr <- predict(modc.PQL3)
pred <- prediction(as.vector(pr), data2$cramps)
auc.perf = performance(pred, measure = "auc")
auc.perf@y.values

pr <- predict(modc.PQL4)
pred <- prediction(as.vector(pr), data2$cramps)
auc.perf = performance(pred, measure = "auc")
auc.perf@y.values

等...

2 つのモデルの精度が異なるかどうかを判断するための基準が設定されていないため (AIC を使用する場合の 2 点ルールなど)、これは適切ではありません。したがって、1 つのモデルの精度値が高くても、すべてのモデルがほぼ同等である可能性があります。

最後に、glmmPQL で半分のノルム プロットを実行する方法を知っている人がいたら教えてください!

多くの感謝を込めて

アレックス

4

0 に答える 0