Eu quero ajustar o modelo misto usando lme4, nlme, pacote de regressão baiano ou qualquer outro disponível.
Modelo misto nas convenções de codificação Asreml-R
Antes de entrar em detalhes, talvez desejemos ter detalhes sobre as convenções asreml-R, para aqueles que não estão familiarizados com os códigos ASREML.
y = Xτ + Zu + e ........................(1) ; o modelo misto usual com, y denota o vetor n × 1 de observações, onde τ é o vetor p × 1 de efeitos fixos, X é uma matriz de projeto n × p de classificação de coluna completa que associa observações à combinação apropriada de efeitos fixos , u é o vetor q × 1 de efeitos aleatórios, Z é a matriz de projeto n × q que associa observações à combinação apropriada de efeitos aleatórios e e é o vetor n × 1 de erros residuais.O modelo (1) é chamado um modelo misto linear ou modelo de efeitos lineares mistos. É assumido
onde as matrizes G e R são funções dos parâmetros γ e φ, respectivamente.
O parâmetro θ é um parâmetro de variação ao qual nos referiremos como o parâmetro de escala.
Nos modelos de efeitos mistos com mais de uma variação residual, surgindo, por exemplo, na análise de dados com mais de uma seção ou variável, o parâmetro θ é fixado em uma. Nos modelos de efeitos mistos com uma única variação residual, então θ é igual à variação residual (σ2). Nesse caso, R deve ser matriz de correlação. Mais detalhes sobre os modelos são fornecidos no manual do Asreml (link) .
Estruturas de variação para os erros: estrutura R e estruturas de variação para efeitos aleatórios: É possível especificar estruturas G.
modelagem de variação em asreml (), é importante entender a formação de estruturas de variação via produtos diretos. A suposição usual de mínimos quadrados (e o padrão em asreml ()) é que eles são distribuídos de forma independente e identificada (IID). No entanto, se os dados fossem de um experimento de campo disposto em uma matriz retangular de r linhas por c colunas, digamos, poderíamos organizar os resíduos e como uma matriz e potencialmente considerar que eles estavam correlacionados automaticamente dentro de linhas e colunas. um vetor em ordem de campo, ou seja, classificando as linhas de resíduos dentro de colunas (plotagens dentro de blocos), a variação dos resíduos pode ser então
são matrizes de correlação para o modelo de linha (ordem r, parâmetro de autocorrelação ½r) e modelo de coluna (ordem c, parâmetro de autocorrelação ½c), respectivamente. Mais especificamente, às vezes é assumida uma estrutura espacial autoregressiva separável bidimensional (AR1 x AR1) para os erros comuns em uma análise de campo.
Os dados de exemplo:
O nin89 é da biblioteca asreml-R, onde diferentes variedades foram cultivadas em réplicas / blocos em campo retangular. Para controlar a variabilidade adicional na direção da linha ou coluna, cada plotagem é referenciada como variáveis de Linha e Coluna (design da coluna da linha). Assim, esta coluna de linha é projetada com bloqueio. O rendimento é medido variável.
Modelos de exemplo
Preciso de algo equivalente aos códigos asreml-R:
A sintaxe do modelo simples terá a seguinte aparência:
rcb.asr <- asreml(yield ∼ Variety, random = ∼ Replicate, data = nin89)
.....model 0O modelo linear é especificado nos argumentos fixo (obrigatório), aleatório (opcional) e rcov (componente de erro) como objetos de fórmula. O padrão é um termo de erro simples e não precisa ser formalmente especificado para o termo de erro, como no modelo 0 .
aqui a variedade é de efeito fixo e a aleatória é replicada (blocos). Além dos termos aleatórios e fixos, podemos especificar o termo do erro. Qual é o padrão neste modelo 0. O componente residual ou de erro do modelo é especificado em um objeto de fórmula através do argumento rcov, consulte os seguintes modelos 1: 4.
O modelo a seguir1 é mais complexo, no qual as estruturas G (aleatória) e R (erro) são especificadas.
Modelo 1:
data(nin89)
# Model 1: RCB analysis with G and R structure
rcb.asr <- asreml(yield ~ Variety, random = ~ idv(Replicate),
rcov = ~ idv(units), data = nin89)Este modelo é equivalente ao modelo 0 acima e introduz o uso do modelo de variação G e R. Aqui, a opção random e rcov especifica as fórmulas random e rcov para especificar explicitamente as estruturas G e R. onde idv () é a função especial do modelo em asreml () que identifica o modelo de variação. A expressão idv (units) define explicitamente a matriz de variação de e para uma identidade em escala.
# Modelo 2: modelo espacial bidimensional com correlação em uma direção
sp.asr <- asreml(yield ~ Variety, rcov = ~ Column:ar1(Row), data = nin89)unidades experimentais de nin89 são indexadas por coluna e linha. Portanto, esperamos variação aleatória em duas direções - direção de linha e coluna neste caso. onde ar1 () é uma função especial que especifica um modelo de variação autoregressiva de primeira ordem para Row. Essa chamada especifica uma estrutura espacial bidimensional para erro, mas com correlação espacial apenas na direção da linha. O modelo de variação para Column é identidade (id ()), mas não precisa ser formalmente especificado, pois esse é o padrão.
# modelo 3: modelo espacial bidimensional, estrutura de erro em ambas as direções
sp.asr <- asreml(yield ~ Variety, rcov = ~ ar1(Column):ar1(Row),
data = nin89)
sp.asr <- asreml(yield ~ Variety, random = ~ units,
rcov = ~ ar1(Column):ar1(Row), data = nin89)semelhante ao modelo 2 acima, no entanto, a correlação é de duas direções - autorregressiva.
Não sei ao certo quanto desses modelos é possível com pacotes R de código aberto. Mesmo que a solução de qualquer um desses modelos seja de grande ajuda. Mesmo que o valor de +50 possa estimular o desenvolvimento desse pacote, será de grande ajuda!
Consulte MAYSaseen forneceu a saída de cada modelo e dados (como resposta) para comparação.
Edições: A sugestão a seguir foi recebida no fórum de discussão de modelos mistos: "Você pode ver os pacotes de covariância regressiva e espacial de David Clifford. O primeiro permite o ajuste de modelos mistos (gaussianos), nos quais é possível especificar a estrutura da matriz de covariância de maneira muito flexível (por exemplo, usei-o para dados de linhagem). O pacote spatialCovariance usa regress para fornecer modelos mais elaborados que o AR1xAR1, mas pode ser aplicável. Você pode ter que corresponder ao autor sobre como aplicá-lo ao seu problema exato. "
lme4. Você pode (a) nos dizer por que precisa fazer isso emlme4vez deasreml-R(b) considerar publicarr-sig-mixed-modelsonde há experiência mais relevante?corStructemnlme(para correlações anisotrópicos) ... Seria bom se você pudesse indicar sucintamente (em palavras ou equações) os modelos estatísticos correspondentes a estas declarações ASREML, já que não estamos todos familiarizados com sintaxe ASREML ...MCMCglmme tenho certeza de que (além de ospatialCovariancemencionado, com o qual não estou familiarizado), a única maneira de fazê-lo em R é definindo newcorStructs - o que é possível, mas não trivial.Respostas:
Você pode ajustar esse modelo com o AD Model Builder. O AD Model Builder é um software gratuito para a construção de modelos não lineares gerais, incluindo modelos de efeitos aleatórios não lineares gerais. Assim, por exemplo, você poderia ajustar um modelo espacial binomial negativo, onde a média e a dispersão excessiva tinham uma estrutura ar (1) x ar (1). Criei o código para este exemplo e o ajustei aos dados. Se alguém estiver interessado, provavelmente é melhor discutir isso na lista em http://admb-project.org
Nota: Existe uma versão R do ADMB, mas os recursos disponíveis no pacote R são um subconjunto do software ADMB independente.
Para este exemplo, é mais fácil criar um arquivo ASCII com os dados, lê-lo no programa ADMB, executar o programa e, em seguida, ler as estimativas de parâmetro etc. de volta ao R para o que você quiser fazer.
Você deve entender que o ADMB não é uma coleção de pacotes, mas uma linguagem para escrever um software não-linear de estimativa de parâmetros. Como eu disse antes, é melhor discutir isso na lista do ADMB, onde todos sabem sobre o software. Depois que estiver pronto e você entender o modelo, poderá postar os resultados aqui. No entanto, aqui está um link para os códigos ML e REML que reuni para os dados do trigo.
http://lists.admb-project.org/pipermail/users/attachments/20111124/448923c8/attachment.zip
fonte
Model 0
ASReml-R
lme4
nlme
fonte
Model 1
ASReml-R
nlme
Veja o truque
fonte
Model 2
ASReml-R
nlme
Trabalhando, ainda não descobri. Pode ser algo assim. Ainda não conseguia descobrir como fazer
rcov=~Column:ar1(Row)comnlmefonte
Model 3
ASReml-R
nlme
Trabalhando, ainda não descobri. Pode ser algo assim. Ainda não conseguia descobrir como fazer
rcov=~ar1(Column):ar1(Row)comnlmeNão consegui descobrir como encaixar os modelos 2 e 3
nlme. Estou trabalhando nisso e atualizarei a resposta quando terminar. Mas incluí a saída dosASReml-Rmodelos 2 e 3 para fins de comparação. Kevin tem uma boa experiência em analisar esses modelos e Ben Bolker tem uma autoridade maravilhosa em modelos mistos. Espero que eles possam nos ajudar nos modelos 2 e 3.fonte