# Maximum Likelihood Ratio with restriction

by user1329307   Last Updated October 10, 2019 02:19 AM

I'm trying to code the likelihood ratio test about two categories for observations from a multinomial experiment with twelve categories; e.g.,

$$H_0: p_1 \leq p2$$ versus $$H_1: p_1>p_2$$.

As a test of concepts, I attempted to compare the power of the likelihood ratio test and a naive use of a binomial proportions test for indepentent samples (ignoring the negative correlation between $$p_1$$ and $$p_2$$).

library(extraDistr)

mnom.ll<-function(x,theta,size,neg=FALSE){
p<-theta
ll<-dmnom(x=x,size=size,prob=p,log=TRUE)
if(sum(p)!=1){
ll<- -Inf
}
ifelse(!neg,ll,-ll)
}

mnom.like.unrestricted<-function(x,theta,size,neg=FALSE){
p<-theta
like<-dmnom(x=x,size=size,prob=p)
if(sum(p)!=1){
like<-0
}
ifelse(!neg,like,-like)
}

mnom.like.restricted<-function(x,theta,size,neg=FALSE){
p<-theta
like<-dmnom(x=x,size=size,prob=p)
if(p[1]!=p[2]|sum(p)!=1){
like<- 0 # penalty
}
ifelse(!neg,like,-like)
}

size=1483
binomial.test<-rep(NA,5000)
lr.test<-rep(NA,5000)
for(i in 1:5000){
survey.dat<-c(rmnom(n=1,size=size,prob=c(0.29,0.26,0.16,0.04,0.03,0.03,0.01,0.02,0.02,0.01,0.01,0.12))) #power
binomial.summary<-binom.test(x = survey.dat[c(1:2)], n = rep(sum(survey.dat),2),alternative = "two.sided")
binomial.test[i]<-binomial.summary\$p.value

mle.unrestricted<-optim(par = survey.dat/size,
fn = mnom.like.unrestricted,
x=survey.dat,
size=size,
neg=TRUE)
mle.restricted<-optim(par =  survey.dat/size,
fn = mnom.like.restricted,
x=survey.dat,
size=size,
neg=TRUE)
mle.restricted

ll.unrestricted<-mnom.ll(x=survey.dat,theta=mle.unrestricted$$par,size = size,neg = FALSE) ll.restricted<-mnom.ll(x=survey.dat,theta=mle.restricted$$par,size = size,neg = FALSE)
lr.test[i]<-1-pchisq(-2*(ll.restricted - ll.unrestricted),df=1)
}
length(which(binomial.test<0.05))/5000
length(which(lr.test<0.05))/5000



The issue is that the penalty in the restricted case isn't working even though it's there and I get the same mles in both the restricted and unrestricted case, which match the estimates based on the sample. The issue is that the restricted mles have a convergence code of 0 and the unrestricted have a convergence code of 1. Is it an issue with my penalty? Is there a quick hack around this issue?

Tags :

## Related Questions

Updated June 03, 2015 22:08 PM

Updated March 04, 2017 19:19 PM

Updated May 18, 2018 16:19 PM

Updated August 27, 2018 09:19 AM

Updated February 16, 2019 01:19 AM