TABLE OF CONTENTS

MetropolisHastings/mh.frailty.spline [ Functions ]

NAME

    mh.frailty.spline --- MH for frailty spline parameters

FUNCTION

Metropolis-Hastings steps for spline parameters for frailty

SYNOPSIS

2673 mh.frailty.spline <- function(hazard, frailty, regression)

INPUTS

    hazard     RCurve for hazard
    frailty    RCurve for frailty
    regression RRegression structure

OUTPUTS

    frailty    Rcurve with updated frailty parameters

SOURCE

2676 {
2677     if(!frailty$hasspline) return(frailty)
2678     sumacc <- 0
2679     nj <- length(frailty$spline.par)
2680     ord2 <- round(frailty$spline.ord) / 2
2681     # base likelihood
2682     baselik <- mklik.spline.frail(frailty$spline.par, hazard, frailty, regression)
2683     # update each parameter separately
2684     for(j in 1:nj){
2685         # get position j and position k which will be used to keep the mean at 1
2686         cand <- frailty$spline.par
2687         k <- j
2688         Ej <- frailty$spline.basisexp[j]
2689         while(j == k | k < ord2 | k > nj - ord2) k <- floor(runif(1, 1, nj + 1))
2690         Ek <- frailty$spline.basisexp[k]
2691         cand[j] <- cand[j] + frailty$spline.tun * rnorm(1, 0, frailty$spline.candsd[j])
2692         newmean <- Ej * (exp(cand[j]) - exp(frailty$spline.par[j]))
2693         newinner <- exp(frailty$spline.par[k]) - newmean / Ek
2694         if(newinner <= 0) next
2695         candk = log(newinner)
2696         cand[k] = candk
2697         # candidate likelihood
2698         candlik <- mklik.spline.frail(cand, hazard, frailty, regression)
2699         thisacc <- acceptreject(baselik, candlik, 1)
2700         if(thisacc){
2701             baselik <- candlik
2702             frailty <- updatespline(frailty, cand)
2703             frailty$spline.fvar <- frailtysplinefvar(frailty)
2704         }
2705         sumacc <- sumacc + thisacc
2706     }
2707     frailty$spline.accept <- sumacc / nj
2708     return(frailty)
2709 }