[R-br] RES: Intervalos de confiança com boot.ci

Pedro Emmanuel Alvarenga Americano do Brasil emmanuel.brasil em gmail.com
Quarta Setembro 23 11:13:03 BRT 2015 

Elias,

Na verdade o w entrou na função apenas porque eu estava tentando reproduzir
o primeiro exemplo da função boot.ci. Eu não tinha intenção de colocar
qualquer peso na função. Mas sem o w a função boot retonrna sempre um erro.
Assim coloquei o w sempre com peso 1, achando que assim daria certo mas
sempre com pesos sem "função". O fato é que agora eu voltei a ter o
problema inicial. A função funciona mas retorna intervalos com valores
iguais a zero. Eu esperava que o intervalo contivesse o valor que retorna a
função originalmente que é

1.000319

e que fosse alguma coisa proxima dos intervalos que retornam a função
original SMR. Assim, os intervalos ainda não fazem sentido para mim.  Eu
ainda não consegui entender o que estou fazendo de errado.

> SMR2 <- function(a,w=1){length(a$Hosp_Death[which(a$Hosp_Death==1)])*w/(length(a$SAPS3Pro2) * mean(a$SAPS3Pro2))*w}> SMR2(a)[1] 1.000319> plot(boo <- boot(a,SMR2,999,stype='w'))> boot.ci(boo)BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 999 bootstrap replicates

CALL : boot.ci(boot.out = boo)

Intervals :
Level      Normal              Basic             Studentized
95%   ( 0,  0 )   ( 0,  0 )   ( 0,  0 )

Level     Percentile            BCa
95%   ( 0,  0 )   ( 0,  0 )
Calculations and Intervals on Original Scale> quantile(boo$t,
c(0.025,0.975))        2.5%        97.5%
0.000000e+00 3.777955e-09

>




Pedro Brasil

Em 23 de setembro de 2015 10:50, Elias Teixeira Krainski <
eliaskrainski em yahoo.com.br> escreveu:

> Considerando a funcao que calcula a estatistica, e' mais apropriado usar
> stype='w' na funcao boot()
>
> O CI fica menos amplo que o obtido por GLM:
>
>
> plot(boo <- boot(a,SMR2,999,stype='w'))
> boot.ci(boo)
>
> m <- glm(Hosp_Death ~ offset(log(SAPS3Pro2)), poisson, data=a)
> (coef.stats <- coef(summary(m)))
>
> hist(boo$t)
> abline(v=exp(coef.stats[1,1] + qnorm(c(0.025,0.5,0.975))*coef.stats[1,2]),
> col=2)
> abline(v=quantile(boo$t, c(0.025,0.975)), col=4)
>
> Elias
>
>
> On 23/09/15 15:40, Pedro Emmanuel Alvarenga Americano do Brasil wrote:
>
> ​ Amigos de R,
>
> Desculpa pela demora, so consegui voltar nessa análise hoje. Bom tentei
> fazer o script recomendado pelo Elias. O pedaço do plot do bootstrap
> funciona, mas o pedaço para gerar o intervalo de confiança retorna um erro.
>
> > library(boot)> SMR2 <- function(a,w=1){length(a$Hosp_Death[which(a$Hosp_Death==1)])*w/(length(a$SAPS3Pro2) * mean(a$SAPS3Pro2))*w}> SMR2(a)[1] 1.000319> plot(boo <- boot(a,SMR2,999,stype='i'))> boot.ci(boo)Error in bca.ci(boot.out, conf, index[1L], L = L, t = t.o, t0 = t0.o,  :
>   estimated adjustment 'w' is infinite
>
>
> Achei muito estranho a função plot do bot funcionar e o boot.ci não
> funcionar. Possivelmente ha necessidade de resolver alguma configuração
> para esse pedação do script.
>
> Pedro Brasil
>
> 2015-09-23 0:19 GMT-03:00 ASANTOS <alexandresantosbr em yahoo.com.br>:
>
>> Mauro, segue o script,
>>
>>
>> dput(a[sample(1:nrow(a),100),c("Hosp_Death","SAPS3Pro2")])
>> structure(list(Hosp_Death = c(0L, 1L, 0L, 1L, 0L, 0L, 0L, 1L,
>> 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
>> 1L, 1L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L,
>> 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L,
>> 1L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L,
>> 1L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 1L, 1L, 0L, 0L, 0L, 0L,
>> 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L), SAPS3Pro2 = c(0.1205,
>> 0.3768, 0.315, 0.3768, 0.0311, 0.0719, 0.1735, 0.75, 0.01, 0.1454,
>> 0.1093, 0.0511, 0.0891, 0.3353, 0.0073, 0.0401, 0.0044, 0.0239,
>> 0.5246, 0.1591, 0.0453, 0.0511, 0.8715, 0.1735, 0.6037, 0.4405,
>> 0.2047, 0.0208, 0.1591, 0.2047, 0.1205, 0.2759, 0.0273, 0.0044,
>> 0.2759, 0.0157, 0.0643, 0.2215, 0.2047, 0.1735, 0.1591, 0.0273,
>> 0.0062, 0.0311, 0.0157, 0.4192, 0.0643, 0.1735, 0.239, 0.4618,
>> 0.0401, 0.0802, 0.0273, 0.0239, 0.0025, 0.0157, 0.0311, 0.0157,
>> 0.0719, 0.239, 0.0574, 0.0157, 0.0085, 0.0085, 0.1093, 0.2215,
>> 0.1591, 0.0719, 0.0044, 0.8563, 0.0453, 0.2047, 0.1454, 0.0273,
>> 0.0891, 0.0802, 0.0239, 0.0135, 0.0802, 0.2571, 0.1735, 0.0157,
>> 0.3768, 0.2759, 0.0891, 0.2047, 0.1735, 0.1454, 0.0643, 0.0311,
>> 0.0453, 0.2759, 0.0239, 0.0181, 0.4405, 0.0085, 0.0311, 0.1205,
>> 0.6579, 0.1093)), .Names = c("Hosp_Death", "SAPS3Pro2"), row.names =
>> c(21875L,
>> 6572L, 48643L, 11869L, 35683L, 48246L, 23919L, 29666L, 3676L,
>> 45549L, 36453L, 14510L, 5082L, 1898L, 41549L, 25481L, 28590L,
>> 38198L, 12822L, 12834L, 33267L, 34088L, 47720L, 30063L, 18326L,
>> 11582L, 11454L, 34960L, 18785L, 11385L, 20605L, 28105L, 25496L,
>> 1607L, 48766L, 36534L, 41868L, 45312L, 37206L, 26927L, 38944L,
>> 21087L, 22343L, 1332L, 11496L, 29485L, 38316L, 4486L, 19757L,
>> 45768L, 33028L, 12205L, 13150L, 41270L, 2780L, 44400L, 19696L,
>> 26015L, 14651L, 39093L, 24905L, 17870L, 35016L, 42851L, 20464L,
>> 33155L, 24924L, 33220L, 15379L, 28989L, 33286L, 34782L, 48534L,
>> 9045L, 20403L, 40222L, 8821L, 31240L, 12465L, 11024L, 24407L,
>> 45729L, 7412L, 22344L, 17737L, 3514L, 2335L, 22491L, 22493L,
>> 14077L, 20346L, 12114L, 36868L, 16431L, 25112L, 31272L, 10877L,
>> 14391L, 8422L, 45130L), class = "data.frame")
>>
>>
>> SMR <- function(obs.var,pred.var){
>>   if(length(obs.var) != length(pred.var)){
>>     stop("Length of pred.var and obs.var differ.")
>>    }
>>    if(any(min(pred.var) <0 | max(pred.var) > 1)){
>>      stop("The individual predicted death must range from 0 to 1.")
>>    }
>>    if(any(levels(as.factor(obs.var)) != c(0,1))){
>>      stop("Observed death variable must be coded as 0 and 1.")
>>    }
>>    O <- length(obs.var[which(obs.var==1)])
>>    E <- length(pred.var) * mean(pred.var)
>>    SMR <- O / E
>>    lowerCL <- O / E * (1 - 1 / (9 * O) - 1.96 / (3 * sqrt(O)))^3 * 100
>>    upperCL <- (O + 1) / E * (1 - (1 / (9 * (O + 1))) + 1.96 / (3 * sqrt(O
>> + 1)))^3 * 100
>>    output <- c(SMR=SMR,lower.Cl=lowerCL,upper.Cl=upperCL)
>>    output
>>  }
>>
>>
>> SMR2 <-
>> function(a,w=1){length(a$Hosp_Death[which(a$Hosp_Death==1)])*w/(length(a$SAPS3Pro2)
>> * mean(a$SAPS3Pro2))*w}
>>
>>
>> boot.ci(boot(a, SMR2, R = 10, stype = "w",parallel = "multicore"))
>>
>> x <- boot.ci(boot(a, SMR2, R = 10, stype = "w",parallel = "snow"))
>> SMR2 <- function(a,w=1)
>>     sum(w*a$Hosp_Death)/sum(w*a$SAPS3Pro2)
>>
>> SMR2(a) ### SMR observada nos dados
>>
>> require(boot)
>>
>> plot(boo <- boot(a,SMR2,999,stype='i'))
>> boot.ci(boo)
>>
>> Em 22/09/2015 23:14, Mauro Sznelwar escreveu:
>>
>>> Alguém poderia passar o código anterior, pois eu apaguei o que chegou
>>> antes.
>>>
>>>
>>>
>>> vc pode fazer assim:
>>>
>>> SMR2 <- function(a,w=1)
>>>       sum(w*a$Hosp_Death)/sum(w*a$SAPS3Pro2)
>>>
>>> SMR2(a) ### SMR observada nos dados
>>>
>>> require(boot)
>>>
>>> plot(boo <- boot(a,SMR2,999,stype='i'))
>>> boot.ci(boo)
>>>
>>> .
>>>
>>>
>>> ---
>>> Este email foi escaneado pelo Avast antivírus.
>>> https://www.avast.com/antivirus
>>>
>>> _______________________________________________
>>> R-br mailing list
>>> R-br em listas.c3sl.ufpr.br
>>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e fornea
>>> cdigo mnimo reproduzvel.
>>>
>>
>> --
>> ======================================================================
>> Alexandre dos Santos
>> Proteção Florestal
>> IFMT - Instituto Federal de Educação, Ciência e Tecnologia de Mato Grosso
>> Campus Cáceres
>> Caixa Postal 244
>> Avenida dos Ramires, s/n
>> Bairro: Distrito Industrial
>> Cáceres - MT                      CEP: 78.200-000
>> Fone: (+55) 65 8132-8112 (TIM)   (+55) 65 9686-6970 (VIVO)
>> e-mails:alexandresantosbr em yahoo.com.br
>>         alexandre.santos em cas.ifmt.edu.br
>> Lattes: http://lattes.cnpq.br/1360403201088680
>> ======================================================================
>>
>> _______________________________________________
>> R-br mailing list
>> R-br em listas.c3sl.ufpr.br
>> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e fornea
>> cdigo mnimo reproduzvel.
>>
>
>
>
> _______________________________________________
> R-br mailing listR-br em listas.c3sl.ufpr.brhttps://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
>
> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forne�a c�digo m�nimo reproduz�vel.
>
>
>
> _______________________________________________
> R-br mailing list
> R-br em listas.c3sl.ufpr.br
> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça
> código mínimo reproduzível.
>
-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20150923/4f5ab637/attachment.html>

Mais detalhes sobre a lista de discussão R-br