<div dir="ltr">Maurício,<div><br></div><div style>Então, se bem entendi o seu problema acho que aí vaí uma solução:</div><div style><br></div><div style><div>n.trat <- 3 # número de tratamentos</div><div>n.rep <- 2 # número de repetições</div>
<div>n.time <- 4 # número de tempos</div><div>n.ind <- n.trat * n.rep # número de indivíduos</div><div>n <- n.trat * n.rep * n.time # número total de registros</div><div>id <- rep(1:n.ind, each = n.time) # identificação do indivíduo </div>
<div>time <- rep(1:n.time, times = n.ind) # identificação do tempo </div><div>trat <- rep(factor(LETTERS[1:n.trat]), each = n.rep*n.time) # tratamentos</div><div><br></div><div>dados <- data.frame(id, time, trat) # estrutura do delineamento</div>
<div><br></div><div>y <- function(x) ifelse(dados$time < 2, rbinom(1, 1, .78), rbinom(1, 1, .35)) ## funcao que simula um ensaio</div><div><br></div><div>nsmc <- 1000 # numero de simulacoes de monte carlo</div><div>
varios <- replicate(nsmc, cbind(dados, y = y(x)), simplify = FALSE) # repete o ensaio nsmc vezes</div><div><br></div><div>modelo <- function( plan ) glm(y ~ trat, data = plan, family = binomial)$coef # funcao que faz o ajuste</div>
<div><br></div><div>ajustes <- do.call(rbind, lapply(varios, modelo)) # varios ajustes</div><div><br></div><div>medias <- apply(ajustes, 2, mean) # médias dos coeficientes</div><div><br></div><div style>Primeiro criei uma função que atribui os resultados ao vetor y, segundo suas premissas. O replicate() repete esse procedimento por nsmc vezes e armazena cada uma das simulações em uma posição da lista.</div>
<div style><br></div><div style>Depois disso é só usar a sua função de ajuste (usei a glm normal) em uma função e aplica com o lapply(). Fiz um artifício de já organizar os coeficientes do modelo numa matriz, uma linha para cada simulação e uma coluna por parâmetro.</div>
<div style><br></div><div style>Espero ter ajudado.</div><div style><br></div><div style>abraço,</div><div style>FH</div></div></div><div class="gmail_extra"><br><br><div class="gmail_quote">2013/4/10 Maurício Lordêlo <span dir="ltr"><<a href="mailto:mslordelo@gmail.com" target="_blank">mslordelo@gmail.com</a>></span><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>Caros,</div><div>No CMR abaixo crio uma matriz de delineamento que estah associada ao nro de tratamentos, nro de repeticoes de cada tratamento e o tempo que cada resposta (neste caso binaria) foi observada.</div>
<div>Suponha que o vetor y seja gerado levando em consideracao que no tempo=1, a probabilidade eh 0.78; nos tempos>1, esta probabilidade eh 0.35.<br></div><div><br></div><div>n.trat=3<span style="white-space:pre-wrap"> </span># número de tratamentos<br>
</div><div>n.rep=2<span style="white-space:pre-wrap"> </span># número de repetições</div><div>n.time=4 <span style="white-space:pre-wrap"> </span># número de tempos</div><div>n.ind=n.trat*n.rep<span style="white-space:pre-wrap"> </span># número de indivíduos</div>
<div>n=n.trat*n.rep*n.time<span style="white-space:pre-wrap"> </span># número total de registros</div><div>id=rep(1:n.ind, each = n.time)<span style="white-space:pre-wrap"> </span># identificação do indivíduo </div><div>
time=rep(1:n.time, times = n.ind)<span style="white-space:pre-wrap"> </span># identificação do tempo </div><div>trat=rep(factor(LETTERS[1:n.trat]), each = n.rep*n.time)<span style="white-space:pre-wrap"> </span># tratamentos</div>
<div><br></div><div>y = matrix(data = NA, nrow = n, ncol = 1, byrow = FALSE,dimnames = NULL)<br></div><div>data1=data.frame(id,time,trat,y)<br></div><div>for(i in 1:n)<br></div><div>{</div><div><span style="white-space:pre-wrap"> </span>y[i]=ifelse(data1[i,2]<2,rbinom(1,1,.78),rbinom(1,1,.35))</div>
<div>}</div><div>data2=data.frame(id,time,trat,y)<br></div><div>data2</div><div><br></div><div>Neste caso tenho formado meu conjunto de dados (data2) já com os valores de y gerados.</div><div>Agora vem a minha dificuldade. Preciso rodar um mesmo modelo repetidas vezes, como por exemplo: </div>
<div><br></div><div>modelo1=geeglm(y~trat, id=id, data=data2, family=binomial) </div><div><br></div><div>e depois obter uma media para cada um dos coeficientes (com seus erros-padrao) obtidos. </div><div>Ou seja, terei outros conjuntos de dados gerados (com o mesmo nro de tratamentos, repeticoes e tempos)</div>
<div>sendo que a unica diferenca estarah nos valores de y. Preciso obter a média dos coeficientes de cada um </div><div>dos modelos rodados a partir destes diferentes y´s.</div><div>Agradeco quem puder auxiliar.</div><div>
Mauricio</div><div><br></div></div>
<br>_______________________________________________<br>
R-br mailing list<br>
<a href="mailto:R-br@listas.c3sl.ufpr.br">R-br@listas.c3sl.ufpr.br</a><br>
<a href="https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br" target="_blank">https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br</a><br>
Leia o guia de postagem (<a href="http://www.leg.ufpr.br/r-br-guia" target="_blank">http://www.leg.ufpr.br/r-br-guia</a>) e forneça código mínimo reproduzível.<br></blockquote></div><br></div>