# R/exemple.R In marqLevAlg: A Parallelized General-Purpose Optimization Based on Marquardt-Levenberg Algorithm

```#' Log-likelihood of a linear mixed model with random intercept
#'
#' @param b numeric vector specifying the parameter's values in the
#' following order : first the fixed effects and then the standard
#' deviation of the random intercept and of the independent error
#' @param Y numeric vector including the dependent outcome vector ordered by subject
#' @param X numeric matrix including the covariates
#' @param ni interger vector giving the number of repeated measures for each subject
#' @return the log-likelihood of the linear mixed model at point b
#' @export

loglikLMM <- function(b,Y,X,ni)
{
beta <- b[1:ncol(X)]
alpha <- b[ncol(X)+1]
sigma <- b[ncol(X)+2]

loglik <- 0
j <- 0

for(i in 1:length(ni))
{
Vi <- alpha^2 * matrix(1,nrow=ni[i],ncol=ni[i]) + sigma^2 * diag(ni[i])

loglik <- loglik - (ni[i]*log(2*pi) + log(det(Vi)) +
t(Y[j+1:ni[i]] - X[j+1:ni[i],]%*%beta) %*% solve(Vi) %*% (Y[j+1:ni[i]] - X[j+1:ni[i],]%*%beta) )/2

j <- j + ni[i]
}

return(as.numeric(loglik))
}

#' Gradient of the log-likelihood of a linear mixed model with random intercept
#'
#' @param b numeric vector specifying the parameter's values in the
#' following order : first the fixed effects and then the standard
#' deviation of the random intercept and of the independent error
#' @param Y numeric vector including the dependent outcome vector ordered by subject
#' @param X numeric matrix including the covariates
#' @param ni interger vector giving the number of repeated measures for each subject
#' @return a vector containing the gradient of the log-likelihood of the linear mixed model at point b
#' @export
{
beta <- b[1:ncol(X)]
alpha <- b[ncol(X)+1]
sigma <- b[ncol(X)+2]

j <- 0

for(i in 1:length(ni))
{
Vi <- alpha^2 * matrix(1,nrow=ni[i],ncol=ni[i]) + sigma^2 * diag(ni[i])
invVi <- solve(Vi)

t(Y[j+1:ni[i]] - X[j+1:ni[i],,drop=FALSE]%*%beta) %*% (-2*alpha* invVi %*% matrix(1,nrow=ni[i],ncol=ni[i]) %*% invVi) %*% (Y[j+1:ni[i]] - X[j+1:ni[i],,drop=FALSE]%*%beta))/2

t(Y[j+1:ni[i]] - X[j+1:ni[i],,drop=FALSE]%*%beta) %*% (-2*sigma* invVi %*% invVi) %*% (Y[j+1:ni[i]] - X[j+1:ni[i],,drop=FALSE]%*%beta))/2

j <- j + ni[i]
}