生態学における一般的な状況は、バイナリ結果 (0 = 死ぬ、1 = 生き残る) を持つ生存モデルであり、個体 (この例では、鳥による個々の営巣の試みを考えてください) は、死亡のリスクにさらされる日数が異なります。これを説明するために、露出日数をリンク関数に組み込む修正ロジスティック回帰を使用します。
Shaffer (2004) の説明:
「毎日の生存率は、適切な予測関数の選択を通じて x に関してモデル化されます。この場合、0 から 1 の間の値が得られるはずです。ロジスティック回帰で行われるように、S 字型のロジスティック関数を使用します。
一般化された線形モデルの体系的なコンポーネントは、[s(x)]t です。次に、関数を検討します。
上記の関数は単調で、θ に関して微分可能であり、g(θ) = β0 + β1x であることが示され、これは一般化線形モデルのリンク関数の基準を満たします。これら 3 つのコンポーネント: 二項応答分布、式 1 で与えられる予測子関数、および式 2 で与えられるリンク関数は、一般化された線形モデルを完全に指定します。このモデル（以下、「ロジスティック暴露モデル」）は、ロジスティック回帰モデルに似ていますが、リンク関数の形式が異なります。ロジスティック露出リンク関数には、ロジスティック回帰リンク関数には存在しない分子と分母の指数 (1/t) が含まれています。指数は、間隔を生き残る確率が間隔の長さに依存するという事実を説明するために必要です。」</p>
このリンク関数のコードは Web で入手できます。また、「help(family)」と入力すると、R で説明されているリンク関数の例の 1 つにもなります。
logexp <- function(days = 1)
{
linkfun <- function(mu) qlogis(mu^(1/days))
linkinv <- function(eta) plogis(eta)^days
mu.eta <- function(eta) days * plogis(eta)^(days-1) *
.Call("logit_mu_eta", eta, PACKAGE = "stats")
valideta <- function(eta) TRUE
link <- paste("logexp(", days, ")", sep="")
structure(list(linkfun = linkfun, linkinv = linkinv,
mu.eta = mu.eta, valideta = valideta, name = link),
class = "link-glm")
}
次のようなモデルで問題なく動作します。
glm(survive ~ date, family=binomial(link=logexp(days=dat$Days)),data=dat)
私が遭遇する問題は、ランダム効果を追加して GLMER モデルでこのカスタム リンク関数を使用しようとしたときです (このメソッドが実装されているオンラインの例を 1 つ見つけました: http://rstudio-pubs-static.s3 .amazonaws.com/4082_51aa699bd9f041c7b3f7cf7b9252f60c.html )。
私たちの場合、サイトをランダム効果として含めたいと思います。モデルは、以前の GLM と同じ方法で定式化されます。
glmer(survive ~ date + (1|site), family=binomial(link=logexp(days=dat$Days)),data=dat)
ただし、エラーメッセージが表示されるようになりました。
famType(glmFit$family) のエラー: 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(4)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(2)' 不明なリンク: 'logexp (3)'不明なリンク:'logexp(3)'不明なリンク:'logexp(4)'不明なリンク:'logexp(3)'不明なリンク:'logexp(2)'不明なリンク:'logexp(1)'不明リンク: 'logexp(4)' 不明なリンク: 'logexp(5)' 不明なリンク: 'logexp(4)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(4)' 不明なリンク: 'logexp( 5) '不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(4)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3)' 不明なリンク: 'logexp(3 )'不明なリンク:'logexp(2)'不明なリンク:'logexp(1)'不明なリンク:'logexp(3)'不明なリンク:'logexp(1)'不明なリンク:'logexp(1)'不明なリンク:'logexp(1)'unknown link: 'logexp(1)'unkn さらに: 警告メッセージ: In if (!(lTyp <- match(family$link, linkNms, nomatch = 0))) stop(gettextf("unknownリンク: %s", : 条件の長さが 1 を超えており、最初の要素のみが使用されます
エラー メッセージには、データのすべての行の不明なリンクが、その巣の訪問 (またはデータの行) の曝露日数に対応する番号とともに一覧表示されます。
例: 最初の 'logexp(3)' は、3 日間の暴露日を持つデータの最初の行に対応します。
このカスタム リンク関数を GLMER モデルで使用できた人はいますか? または、エラーの原因について誰かが考えを持っていない場合は?
＃＃＃＃＃＃アップデート＃＃＃＃＃＃
この問題を解決してくれた Ben Bolker に感謝します。3.0.2 と lme4 の最新バージョンに更新し、リンク関数 Ben の R 関連の投稿 ( https://www.rpubs.com/bbolker/logregexp ) を使用しました。これは次のとおりです。
library(MASS)
logexp <- function(exposure = 1)
{
linkfun <- function(mu) qlogis(mu^(1/exposure))
## FIXME: is there some trick we can play here to allow
## evaluation in the context of the 'data' argument?
linkinv <- function(eta) plogis(eta)^exposure
mu.eta <- function(eta) exposure * plogis(eta)^(exposure-1) *
.Call(stats:::C_logit_mu_eta, eta, PACKAGE = "stats")
valideta <- function(eta) TRUE
link <- paste("logexp(", deparse(substitute(exposure)), ")",
sep="")
structure(list(linkfun = linkfun, linkinv = linkinv,
mu.eta = mu.eta, valideta = valideta,
name = link),
class = "link-glm")
}