estimateLRdistn: Estimate "the" limiting distribution of the likelihood ratio...
In logcondens.mode: Compute MLE of Log-Concave Density on R with Fixed Mode, and Perform Inference for the Mode.

Description Usage Arguments Details Value Author(s) References See Also Examples

Sampling from a given distribution, we estimate via Monte Carlo the limiting distribution of 2-log-likelihood-ratio of the modally-constrained log-concave MLE to the (unconstrained) log-concave MLE.

1 2	estimateLRdistn(rdist = rnorm, mode = 0, N.MC = 1e2, n.SS = 10000, xgrid=NULL, prec = 10^-10, seedVal = NULL, debugging = NULL)

`rdist`	A function taking an integer argument `n` and returning n values simulated from a distribution. The distribution is generally log-concave (otherwise we are in a misspecified setting).
`mode`	fixed/known location of mode for constrained estimator.
`N.MC`	Number of Monte Carlo simulations to do for the limiting distribution.
`n.SS`	Sample Size used for each Monte Carlo. (Each MC simulates `n.SS` values from rdist and computes constrained and unconstrained MLE).
`xgrid`	Governs the generation of weights for observations. `NULL` then data are used as they are. Otherwise can be a single numeric of a numeric vector of length `n.SS`. Please see `preProcess` for details.
`prec`	Precision variable
`seedVal`	An optional seed value
`debugging`	Turns off/on debugging. Any non-character value turns debugging off. If debugging is a character string, then this string gives the name of an output file to which `myxx` (the simulated data from rdist), `myxx.uniq` (the corresponding unique values), and `rdist` are saved. If the code crashes, this can be examined. If debugging is on (i.e., is a character) then if `TLLRs[i]` is less than 0, the value of `myxx` will be saved to a file with name given by `paste(debugging, "tmpxxs", i, ".rsav", sep="")`, along with corresponding weights `myww` and the `mode` passed in.

Computes an estimate of the asymptotic distribution of the likelihood ratio statistic 2 (\mbox{log} \hat{f}_n - \mbox{log} \hat{f}_n^0) under the assumption that the true log-concave density f_0 satisfies f_0''(m)<0 where m is the true mode of f_0. The estimate is computed based on a sample of size n.SS from rdist via N.MC Monte Carlo iterations.

Note: the object LCTLLRdistn was created by output from this function with n.SS set to 1.2e3 and N.MC set to 1e4. Thus, estimateLRdistn is _NOT_ needed to simply compute fairly accurate quantiles of the limit distribution of the likelihood ratio statistic. estimateLRdistn is more useful for research purposes. For instance, by passing to mode values that are not the true mode of myr, the statistic can be studied under the alternative hypothesis.

A list(LRs,TLLRs), i.e., "likelihood ratio" and "two log likelihood ratios". Both are numeric vectors of length N.MC.

Note that theoretically all elements of LRs should be nonnegative, but in practice some rounding errors can occur when n.SS is very large.

Charles Doss cdoss@stat.washington.edu,
http://www.stat.washington.edu/people/cdoss/

Duembgen, L, Huesler, A. and Rufibach, K. (2010) Active set and EM algorithms for log-concave densities based on complete and censored data. Technical report 61, IMSV, Univ. of Bern, available at http://arxiv.org/abs/0707.4643.

Duembgen, L. and Rufibach, K. (2009) Maximum likelihood estimation of a log-concave density and its distribution function: basic properties and uniform consistency. Bernoulli, 15(1), 40–68.

Duembgen, L. and Rufibach, K. (2011) logcondens: Computations Related to Univariate Log-Concave Density Estimation. Journal of Statistical Software, 39(6), 1–28. http://www.jstatsoft.org/v39/i06

Doss, C. R. (2013). Shape-Constrained Inference for Concave-Transformed Densities and their Modes. PhD thesis, Department of Statistics, University of Washington, in preparation.

Doss, C. R. and Wellner, J. A. (2013). Inference for the mode of a log-concave density. Technical Report, University of Washington, in preparation.

See activeSetLogCon and activeSetLogCon.mode, which compute the unconstrained and constrained MLEs, which form the likelihood ratio. The object LCTLLRdistn was created by output from this function.

myseed <- 561

{if(require(distr)){
  mydistn <- Norm() ##demonstrate use of distr package
  myr <- mydistn@r
}
else {
  myr <- rnorm
}}

hypothesis.mode <- 0
N.MC <- 100 ## should increase these values for better estimate
n.SS <- 50

LRres <- estimateLRdistn(rdist=myr, mode=hypothesis.mode, N.MC=N.MC, prec=10^-10,
               n.SS=n.SS, seedVal=myseed,
               debugging=FALSE)
TLLRs <- sort(LRres$TLLRs) ##sort is unnecessary, just for examining data
negIdcs <- TLLRs<=0; ## rounding errors 
Nneg <- sum(negIdcs)
print(Nneg)
TLLRs[negIdcs] <- 0

cdf.empirical.f <- ecdf(TLLRs)
xlims <- c(min(TLLRs), max(TLLRs))
xpts <- seq(from=xlims[1], to=xlims[2], by=.001)
plot(xpts, cdf.empirical.f(xpts), type="l",
     xlab="TLLRs", ylab="Probability")


#### LCTLLRdistn used 1e4 Monte Carlos with 1.2e3 samples each Monte
####Carlo.
##lines(xpts, LCTLLRdistn@p(xpts), col="blue") ## "object
##'C_R_approxfun' not found" error on winbuilder