#Caros,

#Estou fazendo uso da função "bild" pertencente ao pacote com mesmo nome.

#Nos dados que são usados como exemplo do pacote (locust) temos:

install.packages("bild")

require(bild)

locust2 <- bild(move ~ time+feed, data = locust,start = NULL, aggregate = feed, dependence = "MC1R")

summary(locust2)

#O programa retorna os resultados sem problemas, como mostrado abaixo:

Call:

bild(formula = move ~ time + feed, data = locust, aggregate = feed,start = NULL, dependence = "MC1R")

Number of profiles in the dataset:  24 

Number of profiles used in the fit:  24 

Log likelihood:  -1537.923 

AIC:  3085.845 

Coefficients:   

            Label     Value Std. Error t value  p-value

(Intercept)     1 -1.048305 0.30730336  -3.411 0.000647

time            2  1.572124 0.14600334  10.768 0.000000

feed1           3 -3.161945 0.43072021  -7.341 0.000000

log.psi1        4  1.070874 0.09833972  10.890 0.000000

Random effect (omega):  

      Value  Std. Error 

-0.03532051  0.34195846 

Message:  0 

#Porém, se faço

warnings()

#Ele me retorna

Mensagens de aviso:

1: In sqrt(prob * (1 - prob)) : NaNs produzidos

2: In sqrt(prob * (1 - prob)) : NaNs produzidos

3: In sqrt(prob * (1 - prob)) : NaNs produzidos

4: In sqrt(prob * (1 - prob)) : NaNs produzidos

#Mesmo com estas mensagens consigo obter os resultados. 

#O problema surge quando peço para simular muitos ajustes deste modelo.

#Por exemplo,

n.time=4   #numero de medidas no tempo

n.trat=3    #numero de tratamentos

n.rep=10  #numero de repeticoes

ns=15 #numero de simulacoes

n.ind <- n.trat * n.rep  #numero de individuos

n <- n.trat * n.rep * n.time   #numero total de registros

id <- rep(1:n.ind, each = n.time)   #identificacao do individuo

time <- rep(1:n.time, times = n.ind)   #tempo

trat <- rep(factor(LETTERS[1:n.trat]), each = n.rep*n.time) #tratamentos 

dados=data.frame(id,time,trat) #estrutura do banco de dados

#funcao para gerar resposta binaria

z <- function() 

{

 y <- integer(n)

 for(i in 1:n) 

 y[i] <- rbinom(1,1,.5)

 return(y)

}

repet <- replicate(ns, z(), simplify = FALSE)

#funcao para obter as estimativas

coef_bild <- function(y) 

{

dat <- data.frame(dados, y=z())

coefs=bild(y ~ trat+time, data = dat, start=NULL,dependence="MC1R")@coefficients

coefs_bild=coefs[1:4]

}

#lista dos coeficientes ajustados

coefs_ajust <- do.call(rbind, lapply(repet, coef_bild))

#obtencao das medias dos coeficientes estimados

media_coefs <- apply(coefs_ajust, 2, mean)

#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.

#Quando coloco ns=2000 (que é a quantidade que preciso)  demora bastante tempo, não retorna o que preciso e emite a mensagem:

Error in solve.default(Info) :

  system is computationally singular: reciprocal condition number = 1.69179e-16

In addition: There were 50 or more warnings (use warnings() to see the first 50)

#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?

# 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.

#Agradeço antecipadamente,

#Maurício