Aby zrobić ten wykres, wygenerowałem losowe próbki o różnej wielkości z rozkładu normalnego ze średnią = 0 i sd = 1. Przedziały ufności zostały następnie obliczone przy użyciu wartości odcięcia alfa w zakresie od 0,001 do .999 (czerwona linia) za pomocą funkcji t.test (), prawdopodobieństwo profilu zostało obliczone przy użyciu kodu poniżej, który znalazłem w notatkach z wykładów umieszczonych w linii (mogę t znajdź link w tej chwili Edytuj: Znaleziono ), jest to pokazane przez niebieskie linie. Zielone linie pokazują znormalizowaną gęstość za pomocą funkcji gęstości R (), a dane są pokazane przez wykresy prostokątne u dołu każdego wykresu. Po prawej stronie znajduje się wykres gąsienicowy z 95% przedziałami ufności (czerwony) i 1/20 maksymalnych przedziałów prawdopodobieństwa (niebieski).

Kod R użyty dla prawdopodobieństwa profilu:

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

Moje konkretne pytanie dotyczy tego, czy istnieje znany związek między tymi dwoma typami przedziałów i dlaczego przedział ufności wydaje się bardziej zachowawczy dla wszystkich przypadków, z wyjątkiem sytuacji, gdy n = 3. Pożądane są również komentarze / odpowiedzi na temat tego, czy moje obliczenia są prawidłowe (i lepszy sposób na zrobienie tego) oraz ogólny związek między tymi dwoma typami przedziałów.

Kod R:

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)

Nie udzielę pełnej odpowiedzi (trudno mi zrozumieć, co dokładnie robisz), ale postaram się wyjaśnić, w jaki sposób budowane jest prawdopodobieństwo profilu. Mogę uzupełnić odpowiedź później.

Pełne się prawdopodobieństwo prawidłowej próbki o wymiarach jest L ( μ , σ 2 ) = ( σ 2 ) - n / 2 exp ( - Σ I ( x i - μ ) 2 / 2 σ 2 )$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")

Link z prawdopodobieństwem Spróbuję podświetlić link z prawdopodobieństwem na poniższym wykresie.

Najpierw określ prawdopodobieństwo:

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

Następnie wykonaj wykres konturowy:

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)

Wartości prawdopodobieństwa profilu to wartości przyjęte przez prawdopodobieństwo wzdłuż czerwonej paraboli.

$\hat\mu$

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

Możesz również użyć prawdopodobieństwa profilu do zbudowania testów wyników, na przykład.

$\chi^2_k$

$0.147$$95\%$

Są to klasyczne wyniki, dlatego przedstawię kilka odniesień na ten temat:

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

Poniższy kod R pokazuje, że nawet w przypadku małych próbek przedziały uzyskane przy obu podejściach są podobne (ponownie używam przykładu Elvisa):

Pamiętaj, że musisz użyć znormalizowanego prawdopodobieństwa profilu.

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

Jeśli użyjemy większej próby, przedziały ufności będą jeszcze bliższe:

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

WAŻNY PUNKT:

Należy pamiętać, że w przypadku konkretnych próbek różne rodzaje przedziałów ufności mogą się różnić pod względem długości lub lokalizacji, tak naprawdę liczy się ich zasięg. Na dłuższą metę wszystkie powinny zapewniać taki sam zasięg, niezależnie od tego, jak bardzo różnią się dla poszczególnych próbek.

$\chi^2$$normalized$

1. Prawdopodobieństwo dziennika profilu jest w przybliżeniu kwadratowe
2. Istnieje transformacja parametru, która sprawia, że ​​prawdopodobieństwo logarytmu profilu jest w przybliżeniu kwadratowe.

Kwadrat jest ważny, ponieważ określa normalny rozkład w skali logarytmicznej. Im bardziej kwadratowy, tym lepsze przybliżenie i wynikające CI. Twój wybór 1/20 wartości granicznej dla przedziałów prawdopodobieństwa jest równoważny ponad 95% CI w granicy asymptotycznej, dlatego właśnie niebieskie przedziały są na ogół dłuższe niż czerwone.

Jest jeszcze jeden problem z prawdopodobieństwem profilu, który wymaga uwagi. Jeśli masz wiele zmiennych, których profilujesz, to jeśli liczba punktów danych na wymiar jest niska, prawdopodobieństwo profilu może być bardzo stronnicze i optymistyczne. W celu zmniejszenia tego obciążenia wykorzystuje się zatem marginalne, warunkowe i zmodyfikowane prawdopodobieństwa profilu.

Tak więc odpowiedź na twoje pytanie brzmi TAK ... związek jest asymptotyczną normalnością większości estymatorów maksymalnego prawdopodobieństwa, co przejawia się w rozkładzie chi-kwadrat współczynnika prawdopodobieństwa.