# plkhci: Profile-likelihood based confidence intervals In Bhat: General likelihood exploration

## Description

function to find `prob`*100% confidence intervals using profile-likelihood. Numerical solutions are obtained via a modified Newton-Raphson algorithm. The method is described in Venzon and Moolgavkar, Journal of the Royal Statistical Society, Series C vol 37, no.1, 1988, pp. 87-94.

## Usage

 `1` ```plkhci(x, nlogf, label, prob=0.95, eps=.001, nmax=10, nfcn=0) ```

## Arguments

 `x` a list with components 'label' (of mode character), 'est' (the parameter vector with the initial guess), 'low' (vector with lower bounds), and 'upp' (vector with upper bounds) `nlogf` the negative log of the density function (not necessarily normalized) for which samples are to be obtained `label` parameter for which confidence bounds are computed `prob` probability associated with the confidence interval `eps` a numerical value. Convergence results when all (logit-transformed) derivatives are smaller `eps` `nmax` maximum number of Newton-Raphson iterations in each direction `nfcn` number of function calls

## Value

2 component vector giving lower and upper p% confidence bounds

## Note

At this point, only a single parameter label can be passed to plkhci. This function is part of the Bhat exploration tool

## Author(s)

E. Georg Luebeck (FHCRC)

`dfp`, `newton`, `logit.hessian`

## Examples

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24``` ``` # generate some Poisson counts on the fly dose <- c(rep(0,50),rep(1,50),rep(5,50),rep(10,50)) data <- cbind(dose,rpois(200,20*(1+dose*.5*(1-dose*0.05)))) # neg. log-likelihood of Poisson model with 'linear-quadratic' mean: nlogf <- function (x) { ds <- data[, 1] y <- data[, 2] g <- x[1] * (1 + ds * x[2] * (1 - x[3] * ds)) return(sum(g - y * log(g))) } # for example define x <- list(label=c("a","b","c"),est=c(10.,10.,.01),low=c(0,0,0),upp=c(100,20,.1)) # get MLEs using dfp: r <- dfp(x,f=nlogf) x\$est <- r\$est plkhci(x,nlogf,"a") plkhci(x,nlogf,"b") plkhci(x,nlogf,"c") # e.g. 90% confidence bounds for "c" plkhci(x,nlogf,"c",prob=0.9) ```

### Example output

```Sat May 26 16:37:04 2018

starting at
1   fmin:  25053.43     10.00 10.00  0.01

reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
reducing step size
covariance matrix is not positive-definite or nearly singular
dfp search intermediate output:
iter:  0   fmin:  -5722.875    nfcn:  23
start main loop:
DFP search completed with status  0

Optimization Result:
iter:  14   fmin:  -25518.78    nfcn:  178

label estimate    low         high
1 a     19.672695   18.610849   20.779657
2 b     0.52007907  0.45697481  0.59163367
3 c     0.051330289 0.046352266 0.05628206

neg. log. likelihood:  -25518.78

will attempt to compute both bounds (+/- direction)

trying lower bound ------------------------
starting at:    0.9657513
initial guess:  20.22061 0.4935567 0.05059988

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.84169
confidence bound of  a  is  20.77392
log derivatives:     -0.01737574 -0.006629493
label estimate  log deriv   log curv
1 a     20.7739   -55.9705    5563.48
2 b     0.467573  -0.0173757  2883.67
3 c     0.0497386 -0.00662949 454.947

trying upper bound ------------------------
starting at:    0.9557458
initial guess:  19.13606 0.5479867 0.05206013

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.841336
confidence bound of  a  is  18.60534
log derivatives:     -0.019745 -0.003788647
label estimate  log deriv   log curv
1 a     18.6053   55.582      5801.87
2 b     0.576527  -0.019745   3187.37
3 c     0.0526788 -0.00378865 553.572

[1] 18.60534 20.77392
Warning messages:
1: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

2: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

neg. log. likelihood:  -25518.78

will attempt to compute both bounds (+/- direction)

trying lower bound ------------------------
starting at:    0.9830178
initial guess:  19.23707 0.5547335 0.0531601

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.841255
confidence bound of  b  is  0.5900483
log derivatives:     0.001746996 -0.04062658
label estimate log deriv  log curv
1 a     18.793   0.001747   5809.89
2 b     0.590048 -30.365    3205.41
3 c     0.054662 -0.0406266 587.969

trying upper bound ------------------------
starting at:    0.9399685
initial guess:  20.11573 0.4875352 0.0494969

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.841807
confidence bound of  b  is  0.4555413
log derivatives:     0.05348663 -0.06622347
label estimate  log deriv  log curv
1 a     20.5478   0.0534866  5562.12
2 b     0.455541  27.253     2870.3
3 c     0.0472593 -0.0662235 417.611

[1] 0.4555413 0.5900483
Warning messages:
1: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

2: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

neg. log. likelihood:  -25518.78

will attempt to compute both bounds (+/- direction)

trying lower bound ------------------------
starting at:    1.035521
initial guess:  19.51417 0.5453683 0.05381505

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.84129
confidence bound of  c  is  0.0559518
log derivatives:     0.1181582 0.07858222
label estimate  log deriv log curv
1 a     19.3912   0.118158  5724.99
2 b     0.568441  0.0785822 3075.5
3 c     0.0559518 -22.2201  589.781

trying upper bound ------------------------
starting at:    0.8904929
initial guess:  19.83219 0.4959326 0.04883894

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  3.841794
confidence bound of  c  is  0.04592346
log derivatives:     0.2032383 0.1572934
label estimate  log deriv log curv
1 a     20.0279   0.203238  5635.08
2 b     0.46983   0.157293  2983.97
3 c     0.0459235 15.9238   406.846

[1] 0.04592346 0.05595180
Warning messages:
1: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

2: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.

neg. log. likelihood:  -25518.78

will attempt to compute both bounds (+/- direction)

trying lower bound ------------------------
starting at:    0.7204909
initial guess:  19.53959 0.5412233 0.05341652

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  2.705379
confidence bound of  c  is  0.05525405
log derivatives:     0.0874495 0.05734045
label estimate  log deriv log curv
1 a     19.433    0.0874495 5719.05
2 b     0.560745  0.0573404 3070.03
3 c     0.0552541 -18.3131  577.082

trying upper bound ------------------------
starting at:    0.6347411
initial guess:  19.80648 0.4997397 0.04923941

begin Newton-Raphson search for profile lkh conf. bounds:
eps value for stop criterium: 0.001
nmax                        : 10

CONVERGENCE:  5  iterations

chisquare value is:  2.705787
confidence bound of  c  is  0.04685451
log derivatives:     0.1343002 0.1012921
label estimate  log deriv log curv
1 a     19.9672   0.1343    5643.64
2 b     0.477939  0.101292  2993.28
3 c     0.0468545 13.7563   423.373

[1] 0.04685451 0.05525405
Warning messages:
1: In idir * as.vector(dbl0) * sqrt(qchi/(dbb0 - tmp0)) :
Recycling array of length 1 in vector-array arithmetic is deprecated.