<div dir="ltr"><div>#Caros,</div><div>#Estou fazendo uso da função "bild" pertencente ao pacote com mesmo nome.</div><div>#Nos dados que são usados como exemplo do pacote (locust) temos:</div><div>install.packages("bild")</div>
<div>require(bild)</div><div>locust2 <- bild(move ~ time+feed, data = locust,start = NULL, aggregate = feed, dependence = "MC1R")</div><div>summary(locust2)</div><div>#O programa retorna os resultados sem problemas, como mostrado abaixo:</div>
<div>Call:</div><div>bild(formula = move ~ time + feed, data = locust, aggregate = feed,start = NULL, dependence = "MC1R")</div><div>Number of profiles in the dataset: 24 </div><div>Number of profiles used in the fit: 24 </div>
<div>Log likelihood: -1537.923 </div><div>AIC: 3085.845 </div><div>Coefficients: </div><div> Label Value Std. Error t value p-value</div><div>(Intercept) 1 -1.048305 0.30730336 -3.411 0.000647</div>
<div>time 2 1.572124 0.14600334 10.768 0.000000</div><div>feed1 3 -3.161945 0.43072021 -7.341 0.000000</div><div>log.psi1 4 1.070874 0.09833972 10.890 0.000000</div><div>Random effect (omega): </div>
<div> Value Std. Error </div><div>-0.03532051 0.34195846 </div><div>Message: 0 </div><div><br></div><div>#Porém, se faço</div><div>warnings()</div><div>#Ele me retorna</div><div>Mensagens de aviso:</div><div>1: In sqrt(prob * (1 - prob)) : NaNs produzidos</div>
<div>2: In sqrt(prob * (1 - prob)) : NaNs produzidos</div><div>3: In sqrt(prob * (1 - prob)) : NaNs produzidos</div><div>4: In sqrt(prob * (1 - prob)) : NaNs produzidos</div><div>#Mesmo com estas mensagens consigo obter os resultados. </div>
<div>#O problema surge quando peço para simular muitos ajustes deste modelo.</div><div>#Por exemplo,</div><div>n.time=4 #numero de medidas no tempo</div><div>n.trat=3 #numero de tratamentos</div><div>n.rep=10 #numero de repeticoes</div>
<div>ns=15 #numero de simulacoes</div><div>n.ind <- n.trat * n.rep #numero de individuos</div><div>n <- n.trat * n.rep * n.time #numero total de registros</div><div>id <- rep(1:n.ind, each = n.time) #identificacao do individuo</div>
<div>time <- rep(1:n.time, times = n.ind) #tempo</div><div>trat <- rep(factor(LETTERS[1:n.trat]), each = n.rep*n.time) #tratamentos </div><div>dados=data.frame(id,time,trat) #estrutura do banco de dados</div><div>
#funcao para gerar resposta binaria</div><div>z <- function() </div><div>{</div><div> y <- integer(n)</div><div> for(i in 1:n) </div><div> y[i] <- rbinom(1,1,.5)</div><div> return(y)</div><div>}</div><div>repet <- replicate(ns, z(), simplify = FALSE)</div>
<div><br></div><div>#funcao para obter as estimativas</div><div>coef_bild <- function(y) <br></div><div>{</div><div>dat <- data.frame(dados, y=z())</div><div>coefs=bild(y ~ trat+time, data = dat, start=NULL,dependence="MC1R")@coefficients</div>
<div>coefs_bild=coefs[1:4]</div><div>}</div><div>#lista dos coeficientes ajustados</div><div>coefs_ajust <- do.call(rbind, lapply(repet, coef_bild))</div><div>#obtencao das medias dos coeficientes estimados</div><div>media_coefs <- apply(coefs_ajust, 2, mean)</div>
<div><br></div><div><br></div><div>#Com apenas ns=15 ele me dá as mensagens de aviso sem eu pedir (as mesmas do exemplo para o conjunto de dados "locust"), mas me retorna o que preciso.</div><div>#Quando coloco ns=2000 (que é a quantidade que preciso) demora bastante tempo, não retorna o que preciso e emite a mensagem:</div>
<div>Error in solve.default(Info) :</div><div> system is computationally singular: reciprocal condition number = 1.69179e-16</div><div>In addition: There were 50 or more warnings (use warnings() to see the first 50)</div>
<div>#As 50 ou mais mensagens são as mesmas que aparecem quando rodo com uma quantidade menor . Alguém tem alguma dica ou forma de conseguir estes resultados para o número de simulações que desejo?</div><div># Caso nao encontre uma solução, verifiquei que com ns=100 ele faz. Tentei fazer um "for de 1 a 20, com ns=100" das funções que defini como "z" e "coef_bild" e depois concatenar, mas não consegui. Em último caso, peço ajuda a quem saiba fazer para este caso específico.</div>
<div>#Agradeço antecipadamente,</div><div>#Maurício</div></div>