Standardowe algorytmy do wykonywania hierarchicznej regresji liniowej?

Czy istnieją standardowe algorytmy (w przeciwieństwie do programów) do wykonywania hierarchicznej regresji liniowej? Czy ludzie zwykle po prostu wykonują MCMC, czy są bardziej wyspecjalizowane, być może częściowo zamknięte, algorytmy?

Istnieje jeden iteracyjny algorytm uogólnionego najmniejszego kwadratu Harveya Goldsteina (IGLS), a także jego niewielka modyfikacja, ograniczona iteracyjna uogólniona uogólniona najmniejszych kwadratów (RIGLS), który daje obiektywne oszacowania parametrów wariancji.

Algorytmy te są wciąż iteracyjne, więc nie są zamknięte, ale są obliczeniowo prostsze niż MCMC lub maksymalne prawdopodobieństwo. Po prostu iterujesz, aż parametry się zbiegną.

• Goldstein H. Wielopoziomowa mieszana analiza modelu liniowego z wykorzystaniem iteracyjnych uogólnionych najmniejszych kwadratów. Biometrika 1986; 73 (1): 43–56. doi: 10.1093 / biomet / 73.1.43

• Goldstein H. Ograniczał bezstronną iteracyjną uogólnioną estymację metodą najmniejszych kwadratów. Biometrika 1989; 76 (3): 622–623. doi: 10.1093 / biomet / 76.3.622

Aby uzyskać więcej informacji o tym i alternatywach, zobacz np .:

• Stephen W. Raudenbush, Anthony S. Bryk. Hierarchiczne modele liniowe: aplikacje i metody analizy danych . (2. edycja) Sage, 2002.

Pakiet lme4 w R używa iteracyjnie ponownie ważonych najmniejszych kwadratów (IRLS) i karany iteracyjnie ponownie ważonych najmniejszych kwadratów (PIRLS). Zobacz pliki PDF tutaj:

http://rss.acs.unt.edu/Rdoc/library/lme4/doc/index.html

Innym dobrym źródłem „algorytmów obliczeniowych” dla HLM (ponownie w zakresie, w którym postrzegasz je jako podobne specyfikacje jak LMM), byłoby:

• McCulloch, C., Searle, S., Neuhaus, J. (2008). Uogólnione modele liniowe i mieszane. 2. edycja. Wiley. Rozdział 14 - Informatyka.

Wymienione przez nich algorytmy obliczania LMM obejmują:

• Algorytm EM
• Algorytm Newtona Raphsona

Algorytmy, które wymieniają dla GLMM, obejmują:

• Kwadratura numeryczna (kwadratura GH)
• Algorytm EM
• Algorytmy MCMC (jak wspomniałeś)
• Stochastyczne algorytmy aproksymacyjne
• Symulowane maksymalne prawdopodobieństwo

Inne sugerowane przez nich algorytmy GLMM obejmują:

• Karane metody quasi-prawdopodobieństwa
• Przybliżenia Laplace'a
• PQL / Laplace z korekcją błędu początkowego

Jeśli uważasz HLM za rodzaj liniowego modelu mieszanego, możesz rozważyć algorytm EM. Strony 22-23 poniższych notatek kursowych wskazują, jak zaimplementować klasyczny algorytm EM dla modelu mieszanego:

http://www.stat.ucla.edu/~yuille/courses/stat153/emtutorial.pdf

###########################################################
#     Classical EM algorithm for Linear  Mixed Model      #
###########################################################
em.mixed <- function(y, x, z, beta, var0, var1,maxiter=2000,tolerance = 1e-0010)
    {
    time <-proc.time()
    n <- nrow(y)
    q1 <- nrow(z)
    conv <- 1
    L0 <- loglike(y, x, z, beta, var0, var1)
    i<-0
    cat("  Iter.       sigma0                 sigma1        Likelihood",fill=T)
    repeat {
            if(i>maxiter) {conv<-0
                    break}
    V <- c(var1) * z %*% t(z) + c(var0) * diag(n)
    Vinv <- solve(V)
    xb <- x %*% beta
    resid <- (y-xb)
    temp1 <- Vinv %*% resid
    s0 <- c(var0)^2 * t(temp1)%*%temp1 + c(var0) * n - c(var0)^2 * tr(Vinv)
    s1 <- c(var1)^2 * t(temp1)%*%z%*%t(z)%*%temp1+ c(var1)*q1 -
                                                c(var1)^2 *tr(t(z)%*%Vinv%*%z)
    w <- xb + c(var0) * temp1
    var0 <- s0/n
    var1 <- s1/q1
    beta <- ginverse( t(x) %*% x) %*% t(x)%*% w
    L1 <- loglike(y, x, z, beta, var0, var1)
    if(L1 < L0) { print("log-likelihood must increase, llikel <llikeO, break.")
                             conv <- 0
break
}
    i <- i + 1
    cat("  ", i,"  ",var0,"  ",var1,"  ",L1,fill=T)
    if(abs(L1 - L0) < tolerance) {break}  #check for convergence
    L0 <- L1
    }
list(beta=beta, var0=var0,var1=var1,Loglikelihood=L0)
}

#########################################################
#  loglike calculates the LogLikelihood for Mixed Model #
#########################################################
loglike<- function(y, x, z, beta, var0, var1)
}
{
n<- nrow(y)
V <- c(var1) * z %*% t(z) + c(var0) * diag(n)
Vinv <- ginverse(V)
xb <- x %*% beta
resid <- (y-xb)
temp1 <- Vinv %*% resid
(-.5)*( log(det(V)) + t(resid) %*% temp1 )
}