Um dieses Diagramm zu erstellen, habe ich Zufallsstichproben unterschiedlicher Größe aus einer Normalverteilung mit Mittelwert = 0 und sd = 1 generiert. Die Konfidenzintervalle wurden dann unter Verwendung von Alpha-Grenzwerten zwischen 0,001 und 0,999 (rote Linie) mit der Funktion t.test () berechnet. Die Profilwahrscheinlichkeit wurde unter Verwendung des Codes berechnet, unter dem ich in den online gestellten Vorlesungsunterlagen gefunden habe (ich kann ' t den Link im Moment finden Edit: es gefunden ), wird dies durch die blauen Linien dargestellt. Grüne Linien zeigen die normalisierte Dichte unter Verwendung der Funktion R density () und die Daten werden durch die Boxplots am unteren Rand jedes Diagramms angezeigt. Auf der rechten Seite ist eine Raupendiagramm der 95% -Konfidenzintervalle (rot) und 1/20 der maximalen Wahrscheinlichkeitsintervalle (blau) dargestellt.

R Code für die Profilwahrscheinlichkeit:

  #mn=mean(dat)
  muVals <- seq(low,high, length = 1000)
  likVals <- sapply(muVals,
                    function(mu){
                      (sum((dat - mu)^2) /
                         sum((dat - mn)^2)) ^ (-n/2)
                    }
  )

Meine spezifische Frage ist, ob es eine bekannte Beziehung zwischen diesen beiden Intervalltypen gibt und warum das Konfidenzintervall für alle Fälle konservativer erscheint, außer wenn n = 3 ist. Kommentare / Antworten darüber, ob meine Berechnungen gültig sind (und eine bessere Möglichkeit, dies zu tun), und die allgemeine Beziehung zwischen diesen beiden Intervalltypen sind ebenfalls erwünscht.

R-Code:

samp.size=c(3,4,5,10,20,1000)
cnt2<-1
ints=matrix(nrow=length(samp.size),ncol=4)
layout(matrix(c(1,2,7,3,4,7,5,6,7),nrow=3,ncol=3, byrow=T))
par(mar=c(5.1,4.1,4.1,4.1))
for(j in samp.size){


  #set.seed(200)
  dat<-rnorm(j,0,1)
  vals<-seq(.001,.999, by=.001)
  cis<-matrix(nrow=length(vals),ncol=3)
  cnt<-1
  for(ci in vals){
    x<-t.test(dat,conf.level=ci)$conf.int[1:2]
    cis[cnt,]<-cbind(ci,x[1],x[2])
    cnt<-cnt+1
  }


  mn=mean(dat)
  n=length(dat)
  high<-max(c(dat,cis[970,3]), na.rm=T)
  low<-min(c(dat,cis[970,2]), na.rm=T)
  #high<-max(abs(c(dat,cis[970,2],cis[970,3])), na.rm=T)
  #low<--high


  muVals <- seq(low,high, length = 1000)
  likVals <- sapply(muVals,
                    function(mu){
                      (sum((dat - mu)^2) /
                         sum((dat - mn)^2)) ^ (-n/2)
                    }
  )


  plot(muVals, likVals, type = "l", lwd=3, col="Blue", xlim=c(low,high),
       ylim=c(-.1,1), ylab="Likelihood/Alpha", xlab="Values",
       main=c(paste("n=",n), 
              "True Mean=0 True sd=1", 
              paste("Sample Mean=", round(mn,2), "Sample sd=", round(sd(dat),2)))
  )
  axis(side=4,at=seq(0,1,length=6),
       labels=round(seq(0,max(density(dat)$y),length=6),2))
  mtext(4, text="Density", line=2.2,cex=.8)

  lines(density(dat)$x,density(dat)$y/max(density(dat)$y), lwd=2, col="Green")
  lines(range(muVals[likVals>1/20]), c(1/20,1/20), col="Blue", lwd=4)
  lines(cis[,2],1-cis[,1], lwd=3, col="Red")
  lines(cis[,3],1-cis[,1], lwd=3, col="Red")
  lines(cis[which(round(cis[,1],3)==.95),2:3],rep(.05,2), 
        lty=3, lwd=4, col="Red")
  abline(v=mn, lty=2, lwd=2)
  #abline(h=.05, lty=3, lwd=4, col="Red")
  abline(h=0, lty=1, lwd=3)
  abline(v=0, lty=3, lwd=1)

  boxplot(dat,at=-.1,add=T, horizontal=T, boxwex=.1, col="Green")
  stripchart(dat,at=-.1,add=T, pch=16, cex=1.1)

  legend("topleft", legend=c("Likelihood"," Confidence Interval", "Sample Density"),
         col=c("Blue","Red", "Green"), lwd=3,bty="n")

  ints[cnt2,]<-cbind(range(muVals[likVals>1/20])[1],range(muVals[likVals>1/20])[2],
                     cis[which(round(cis[,1],3)==.95),2],cis[which(round(cis[,1],3)==.95),3])
  cnt2<-cnt2+1
}
par(mar=c(5.1,4.1,4.1,2.1))


plot(0,0, type="n", ylim=c(1,nrow(ints)+.5), xlim=c(min(ints),max(ints)), 
     yaxt="n", ylab="Sample Size", xlab="Values")
for(i in 1:nrow(ints)){
  segments(ints[i,1],i+.2,ints[i,2],i+.2, lwd=3, col="Blue")
  segments(ints[i,3],i+.3,ints[i,4],i+.3, lwd=3, col="Red")
}
axis(side=2, at=seq(1.25,nrow(ints)+.25,by=1), samp.size)

Ich werde keine vollständige Antwort geben (es fällt mir schwer zu verstehen, was Sie genau tun), aber ich werde versuchen zu klären, wie die Profilwahrscheinlichkeit aufgebaut ist. Ich kann meine Antwort später vervollständigen.

Die volle Wahrscheinlichkeit für eine normale Stichprobe der Größe ist L ( μ , σ 2 ) = ( σ 2 ) - n / 2 exp ( - ∑ i ( x $n$

$\mu$$\sigma^2$$\mu$

$\exp(-n/2)$

> data(sleep)
> difference <- sleep$extra[11:20]-sleep$extra[1:10]
> Lp <- function(mu, x) {n <- length(x); mean( (x-mu)**2 )**(-n/2) }
> mu <- seq(0,3, length=501)
> plot(mu, sapply(mu, Lp, x = difference), type="l")

Verknüpfung mit der Wahrscheinlichkeit Ich werde versuchen, die Verknüpfung mit der Wahrscheinlichkeit mit der folgenden Grafik hervorzuheben.

Definieren Sie zuerst die Wahrscheinlichkeit:

L <- function(mu,s2,x) {n <- length(x); s2**(-n/2)*exp( -sum((x-mu)**2)/2/s2 )}

Dann mache ein Konturdiagramm:

sigma <- seq(0.5,4, length=501)
mu <- seq(0,3, length=501)

z <- matrix( nrow=length(mu), ncol=length(sigma))
for(i in 1:length(mu))
  for(j in 1:length(sigma))
    z[i,j] <- L(mu[i], sigma[j], difference)

# shorter version
# z <- outer(mu, sigma, Vectorize(function(a,b) L(a,b,difference)))

contour(mu, sigma, z, levels=c(1e-10,1e-6,2e-5,1e-4,2e-4,4e-4,6e-4,8e-4,1e-3,1.2e-3,1.4e-3))

$\widehat{\sigma^2}(\mu)$

hats2mu <- sapply(mu, function(mu0) mean( (difference-mu0)**2 ))
lines(mu, hats2mu, col="red", lwd=2)

Die Werte der Profilwahrscheinlichkeit sind die Werte, die von der Wahrscheinlichkeit entlang der roten Parabel genommen werden.

$\hat\mu$

Für Ihr Konfidenzintervall weichen die Ergebnisse aufgrund der Krümmung der Funktion geringfügig ab$\widehat{\sigma^2}(\mu)$ , aber solange Sie sich nur mit einem kurzen Segment davon befassen, ist es fast linear und der Unterschied ist sehr gering .

Sie können die Profilwahrscheinlichkeit beispielsweise auch zum Erstellen von Bewertungstests verwenden.

$\chi^2_k$

$0.147$$95\%$

Dies sind klassische Ergebnisse und deshalb werde ich hier nur einige Hinweise geben:

http://www.jstor.org/stable/2347496

http://www.stata-journal.com/sjpdf.html?articlenum=st0132

http://www.unc.edu/courses/2010fall/ecol/563/001/docs/lectures/lecture11.htm

http://en.wikipedia.org/wiki/Likelihood-ratio_test

http://en.wikipedia.org/wiki/Likelihood_function#Profile_likelihood

Der folgende R-Code zeigt, dass die Intervalle, die mit beiden Ansätzen erhalten werden, selbst für kleine Stichproben ähnlich sind (ich verwende das Elvis-Beispiel erneut):

Beachten Sie, dass Sie die normalisierte Profilwahrscheinlichkeit verwenden müssen.

data(sleep)
x <- sleep$extra[11:20]-sleep$extra[1:10]
n <- length(x)
Rp <- function(mu) {mean( (x-mean(x))^2 )^(n/2)/mean( (x-mu)^2 )^(n/2) }
Rp(mean(x))

mu <- seq(0,3, length=501)
plot(mu, sapply(mu, Rp), type="l")


Rpt<- function(mu) Rp(mu)-0.147 # Just an instrumental function

# Likelihood-confidence interval of 95% level

c(uniroot(Rpt,c(0.5,1.5))$root,uniroot(Rpt,c(1.51,3))$root)

# t confidence interval

t.test(x,conf.level=0.95)$conf.int

Wenn wir eine größere Stichprobe verwenden, sind die Konfidenzintervalle noch enger:

set.seed(123)
x <- rnorm(100)
n <- length(x)
Rp <- function(mu) {mean( (x-mean(x))^2 )^(n/2)/mean( (x-mu)^2 )^(n/2) }
Rp(mean(x))

mu <- seq(-0.5,0.5, length=501)
plot(mu, sapply(mu, Rp), type="l")


Rpt<- function(mu) Rp(mu)-0.147 # Just an instrumental function

# Likelihood-confidence interval of 95% level

c(uniroot(Rpt,c(-0.4,0))$root,uniroot(Rpt,c(0,0.4))$root)

# t confidence interval

t.test(x,conf.level=0.95)$conf.int

Ein wichtiger Punkt:

Beachten Sie, dass für bestimmte Stichproben unterschiedliche Arten von Konfidenzintervallen in Bezug auf Länge oder Position unterschiedlich sein können. Entscheidend ist jedoch die Abdeckung. Auf lange Sicht sollten alle die gleiche Abdeckung bieten, unabhängig davon, wie stark sie sich für bestimmte Proben unterscheiden.

$\chi^2$$normalized$

1. Die Profil-Log-Wahrscheinlichkeit ist ungefähr quadratisch
2. Es gibt eine Parametertransformation, die die Profilprotokollwahrscheinlichkeit annähernd quadratisch macht.

Das Quadrat ist wichtig, weil es eine Normalverteilung in logarithmischer Skala definiert. Je quadratischer es ist, desto besser sind die Approximation und die daraus resultierenden CIs. Die Auswahl des 1/20-Cutoff für die Wahrscheinlichkeitsintervalle entspricht einem CI von mehr als 95% in der asymptotischen Grenze, weshalb die blauen Intervalle im Allgemeinen länger sind als die roten.

Nun gibt es ein weiteres Problem mit der Profilwahrscheinlichkeit, das einige Aufmerksamkeit erfordert. Wenn Sie viele Variablen haben, über die Sie ein Profil erstellen, und die Anzahl der Datenpunkte pro Dimension gering ist, kann die Profilwahrscheinlichkeit sehr voreingenommen und optimistisch sein. Grenz-, bedingte und modifizierte Profilwahrscheinlichkeiten werden dann verwendet, um diese Verzerrung zu verringern.

Die Antwort auf Ihre Frage lautet also JA ... der Zusammenhang ist die asymptotische Normalität der meisten Schätzer für die maximale Wahrscheinlichkeit, die sich in der Chi-Quadrat-Verteilung des Wahrscheinlichkeitsverhältnisses manifestiert.