Modelo beta Weibull modificada em Análise de
Sobrevivência
Giovana O. Silva∗
Departamento de Estatı́stica - UFBA
Valdemiro Piedade Vigas
Bolsista de Iniciação Cientifı́ca - Departamento de Estatı́stica - UFBA
Resumo
Em aplicações na área de análise de sobrevivência é frequente a ocorrência de função de taxa de
falha em forma de U. As distribuições tradicionalmente usadas não apresentam essa propriedade. Considerando dados censurados, é proposto utilizar a distribuição beta Weibull modificada, que apresenta
função de taxa de falha em forma de U, para descrever o comportamento do tempo de sobrevivência e é
apropriada na discriminação entre alguns modelos probabilı́sticos alternativos, tais como, a distribuição
Weibull, a distribuição Weibull exponenciada, a distribuição Rayleigh generalizada, distribuição Weibull
modificada e a distribuição Weibull modificada generalizada. O método de máxima verossimilhança
será usado para estimar os parâmetros dos modelos. Para selecionar um modelo que melhor se ajusta
aos dados será usado o teste da razão de verossimilhança e os critérios AIC e BIC. Conjunto de dados
real será usado para ilustrar a aplicação da nova distribuição. O software R foi usado para análise
dos dados. Para ilustrar o uso desta distribuição, foi usado conjunto de dados com 148 crianças contaminadas pelo HIV por via vertical do Hospital das Clı́nicas da Faculdade de Medicina de Ribeirão
Preto.
Palavras-chaves: Distribuição beta Weibull modificada; censura; máxima verossimilhança.
1
Introdução
Considere um conjunto de dados cuja variável resposta é o tempo até a ocorrência de um evento de
interesse. Esse tempo é denotado por tempo de sobrevivência ou tempo de falha e o conjunto de observações por dados de sobrevivência. Uma caracterı́stica destes dados é a censura, que é a observação
parcial do tempo de interesse. Existem diversos tipos de censura, neste trabalho considerou-se censura à
direita e mecanismo aleatório, além de que a censura é não informativa. A área da estatı́stica que engloba
metodologias para esse tipo de dados é denominada análise de sobrevivência. Para mais detalhes ver, por
exemplo, Colosimo e Giolo (2006)
As distribuições de probabilidade tradicionalmente usadas em análise de sobrevivência são: exponencial, Weibull, log-normal e log-logı́stica. Entretanto, em diversas aplicações na área de sobrevivência o
tempo de sobrevivência apresenta função de taxa de falha na forma de U e os modelos de probabilidade
citados anteriormente não apresentam este tipo de função. Diante deste fato, o objetivo desta trabalho
é propor o uso da distribuição de probabilidade beta Weibull modificada, proposta por Silva, Ortega e
Cordeiro (2010), para descrever o comportamento do tempo de sobrevivência.
1
e-mail: [email protected]
1
Esta distribuição apresenta flexibilidade para acomodar várias formas da função de taxa de falha, por
exemplo, forma de U , e é apropriada na discriminação entre alguns modelos probabilı́sticos alternativos,
tais como, a distribuição Weibull, a distribuição Weibull exponenciada, a distribuição Rayleigh generalizada, distribuição Weibull modificada e a distribuição Weibull modificada generalizada apresentada por
Carrasco, Ortega e Cordeiro (2008). Assim, espera-se que esta distribuição seja útil em uma variedade de
aplicações nas áreas de confiabilidade, medicina e biologia, bem como em outras áreas de investigação.
Para encontrar as estimativas dos parâmetros dessa distribuição utilizou-se o método de máxima
verossimilhança em que o método numérico do tipo quase Newton foi necessário. Os critérios de seleção
de modelos: Critério de informação de Akaike (AIC) e Critério de informação de Bayes (BIC)e o teste
da razão de verossimilhanças foram usados para selecionar um modelo de probabilidade que melhor se
ajusta aos dados. O sof tware R foi usado para análise dos dados. Para ilustrar o uso desta distribuição,
foi usado conjunto de dados com 148 crianças contaminadas pelo HIV por via vertical do Hospital das
Clı́nicas da Faculdade de Medicina de Ribeirão Preto (Silva, 2004).
2
Distribuição beta Weibull modificada
A distribuição beta Weibull modicada (BWM) é composta por 5 parâmetros e apresenta algumas caracterı́sticas importantes como flexibilidade para acomodar várias formas da função de taxa de falha ( por
exemplo forma de U e unimodal ), é apropriada na discriminação entre alguns modelos probabilisticos
alternativos como a Weibull, Weibull modificada e a exponencial de acordo com os valores dos seus respectivos parâmetros. A função densidade de probabilidade, a função de sobrevivência e função de taxa
de falha são descritas por:
f (t) =
αxγ−1 (γ + λx) exp(λx)
[1 − exp{−αtγ exp(λt)}]α−1 exp{−bαtγ exp(λt)}, t > 0
B(a, b)
S(t) = [1 − I1−exp(−αxγ exp(λx)) (a, b)],
αxγ−1 (γ + λx) exp(λx)
[1 − exp{−αtγ exp(λt)}]α−1 exp{−bαtγ exp(λt)}, t > 0
B(a, b)[1 − I1−exp{−αtγ exp(λt)} (a, b)]
Ry
em que Iy (a, b) = By (a, b)/B(a, b) é a função razão beta incompleta, By (a, b) = 0 wa−1 (1 − w)b−1 dw é a
função beta incompleta e B(a, b) é a função beta.
A distribuição BWM contém como caso especiais várias distribuições conhecidas. Por exemplo, quando
λ = 0 simplica para a distribuição beta Weibull. Se γ = 1 além de λ = 0, obtém-se a distribuição beta
exponencial. A distribuição Weibull modificada generalizada é um caso especial quando b = 1. Se a = 1
além de b = 1, a distribuição Weibull modificada é obtida. Para b = 1 e λ = 0, a distribuição BWM
reduz-se a distribuição Weibull exponenciada. Se γ = 1 além de b = 1 e λ = 0, a distribuição BWM
torna-se a distribuição exponencial exponenciada. Para γ = 2, λ = 0 e b = 1, a distribuição BWM
torna-se distribuição Rayleigh generalizada. A distribuição Weibull é um caso especial quando a = b = 1
e λ = 0.
λ(t) =
3
Estimador de máxima verossimilhança
2
Existem vários métodos que podem ser utilizados para estimar os parâmetros dos modelos probabilı́sticos, neste trabalho foi usado o método de máxima verossimilhança. Considerando que os tempos de
sobrevivência e de censura são independentes e que a censura é não informativa, a função de verossimilhança para todos os mecanismos de censura à direita é dada por:
Y
L(θ)α
[f (xi ; θ)]δi [S(xi ; θ)]1−δi
i=1
em que δ é a variável indicadora de censura e θ é o vetor de parâmetros.
Seja x1 , . . . , xn uma amostra aleatória ¡de ¢tamanho n da distribuição BWM(a, b, α, γ, λ). Também
γ
considere a seguinte reparametrização α = α11 . Então, o logaritmo da função de verossimilhança para
o vetor de parâmetros θ = (a, b, α1 , γ, λ)T é escrito como
o
£
¤ X n£
¤
l(θ) = −r log B(a, b) +
log(γ + λxi ) + log(vi ) − log(xi )
+(a − 1)
X
i∈F
n
o
X
X
log 1 − exp(−vi ) − b
vi +
log(1 − I[1−(exp(−vi ))] (a, b)),
i∈F
¡
¢γ
i∈F
(1)
i∈C
xγi
em que vi = α11
exp(λxi ) e r é o número de observações não censuradas.
O estimador de máxima verossimilhança para o vetor de parâmetros é calculado maximizando a
função de verossimilhança ou, equivalentemente, maximizando o logaritmo da função de verossimilhança.
Entretanto, neste caso, o sistema de equações é não linear, tornado-se necessário usar um algoritmo de
otimização, por exemplo, quase-Newton. Neste trabalho, as estimativas foram obtidas diretamente usando
o sof tware R (comando optim).
No caso em que o tamanho da amostra é grande e sob certas condições de regularidade para a função
de verossimilhança, intervalos de confiança e testes de hipóteses podem ser obtidos usando o fato de que a
√ b
− θ) ∼ Np (0, I(θ)−1 ),
distribuição assintótica para os estimadores de máxima verossimilhança é n(θ
em que I(θ) é a matriz de informação de Fisher e p é o número de parâmetros do modelo (Sen e Singer,
1993).
Como é complicado calcular
a matriz
de informação de Fisher devido às observações censuradas, pode½
¾
∂ 2 l(θ )
b que é um estimador consistente para matriz de
se usar a matriz L̈(θ) = −
, avaliada em θ = θ,
T
∂ θ∂ θ
b é
informção de Fisher (Mudholkar, Srivastava e Kollia, 1996). Assim, a distribuição assintótica para θ
√ b
−1
dada por n(θ − θ) ∼ Np {0; L̈(θ) }, em que p é o número de parâmetros do modelo.
4
Seleção de modelos
Para decidir qual modelo é mais adequado para modelar o tempo de sobrevivência, os critérios de seleção
AIC e BIC foram usados, são dados por: AIC = −2log(L) + 2p e BIC = −2log(L) + plog(n), em que
L representa função de verossimilhança, p representa o número de parâmetros e n representa o número
de observações. Tanto o AIC como o BIC selecionam o modelo que apresenta o menor valor.
4.1
Teste da razão de verossimilhanças
Outra forma de selecionar um modelo é usar o teste da razão de verossimilhanças. A hipóteses testadas
são: H0 : modelo de interesse é adequado versus H1 : modelo de interesse não é adequado. A estatı́stica
de teste é escrita da seguinte forma :
b1 ) − log(L)(θ
b0 )],
T RV = 2[log(L)(θ
3
b0 é o estimador de máxima verossimilhança para modelo de interesse e θ
b1 é o estimador de máxima
em que θ
verossimilhança para modelo mais geral. Esta estatı́stica sob H0 tem aproximadamente distribuiçao quiquadrado e graus de liberdade fornecido pela diferença do número de parâmetros dos modelos comparados.
5
Aplicação
Para ilustrar a aplicação da distribuição beta Weibull modificada foi usado o conjunto de dados que está
relacionado ao tempo (dias) até a ocorrência de sororeversão de 148 crianças contaminadas ao HIV por via
vertical no Hospital das Clı́nicas da Faculdade de medicina de Ribeirão Preto, de 1986 a 2001, sendo que
56,09% dessas observações correspondem à censura. Mais detalhe sobre esse conjunto de dados pode ser
encontrado em Silva, 2004.Vale explicar que sororeversão é o processo de desaparecimento dos anticorpos
anti-HIV do sangue (negativação da sorologia anti-HIV) em um indivı́duo que anteriormente apresentava
sorologia anti-HIV positivo.
O objetivo é saber qual o modelo paramétrico que se apresenta mais apropriado para descrever os
dados. No primeiro momento, foram encontradas as estimativas dos parâmetros da distribuição BWM e
alguns de seus sub-modelos utilizando o método de máxima verossimilhança por meio do comando optim
do Sof tware R.
Essas estimativas com os correspondentes erros padrão, entre parênteses, são dados na Tabela 1 para
diversos modelos ajustados.
Tabela 1: Resultados das estimativas e do erro padrão (entr parênteses) da distribuição beta Weibull
modificada e alguns de seus sub-modelos para o conjunto de dados de sororeversão.
Distribuições
exponencial
Weibull
Weibull modificada
Weibull modificada generalizada
beta Weibull modificada
α
372,1094
177,8135
307,6221
12,54027
307,6217
11,9951
364,2627
16,8264
307,4936
0.0024
γ
-
λ
-
a
-
b
-
3,1114
0.0087
3,1125
0,0086
7,8409
0,0000
7,9685
0,0000
-
-
-
0,0003
0,0000
0.0004
0,0000
0,0051
0,0000
-
-
0,3132
0,0002
0,1937
0,0000
0,0851
0,0000
Para identificar qual dessas distribuições é mais apropriada para modelar o tempo de sobrevivência
desse conjunto de dados, foram usados gráficos que relacionam as estimativas das funções de sobrevivência
das distribuições Weibull, Weibull modificada, Weibull modificada generalizada e beta Weibull modificada
com as estimativas da função de sobrevivência obtidas pelo método de Kaplan-Meier.
Dando prosseguimento, foram aplicados os critérios de seleção de modelo AIC e BIC e o teste da razão
de verossimilhança. A tabela 2 apresenta os resultados referentes a estes métodos.
Da Tabela 2, observa-se que das diferentes distribuições ajustadas aos dados, a distribuição beta
Weibull modificada apresenta um melhor ajuste comparado aos outros modelos probabilı́sticos alternativos, pois o critério AIC e BIC apresentam menores valores em relaçao as outras distribuições consideradas.
4
0.4
0.8
S(t):Kaplan−Meier
0.0
0.4
0.8
1.0
0.6
0.2
0.0
0.0
S(t):Kaplan−Meier
0.4
S(t):Beta Weibul Modificada
0.8
1.0
0.6
0.0
0.2
0.4
S(t):Weibul Modificada Generalizada
0.6
0.0
0.2
0.4
S(t):Weibul Modificada
0.8
1.0
0.8
1.0
0.8
0.6
S(t):Weibul
0.4
0.2
0.0
0.0
0.4
0.8
S(t):Kaplan−Meier
0.0
0.4
0.8
S(t):Kaplan−Meier
Figura 1: Função de sobrevivência estimada pelo método de Kaplan-Meier versus a função de sobrevivência estimada para as seguintes distribuições: Weibull, Weibull modificada, Weibull modificada generalizada e beta Weibull modificada dos dados de sororeversão
Tabela 2: AIC, BIC e logaritmo da função de verossimilhança para o modelo beta Weibull modificada e
alguns dos seus submodelos para o conjunto de dados de sororeversão
Distribuições
Weibull
Weibull modificada
Weibull modificada generalizada
beta Weibull modificada
5
AIC
807,9872
809,6044
795,0154
780,6156
BIC
813,9816
818,5960
807,0042
795,6017
e
l(θ)
401,9936
401,8022
393,5077
385,3078
Pode-se ainda usar o teste razão de verossimilhança para testar se o modelo beta Weibull modificada
ajusta bem os dados. Considere as seguintes hipóteses: H0 : b = 1 versus H1 : H0 não é verdade, isto é,
compara o modelo Weibull modificada generalizada com o modelo BWM. Para este teste de hipóteses, o
valor da estatı́stica de teste (T RV =0,0424) (p-valor < 0.001) que fornece indicações favoráveis ao modelo
BWM.
6
Conclusão
O modelo probabilı́stico com cinco parâmetros, denominado de distribuição beta Weibull modificada, que
estende várias distribuições amplamente abordadas na literatura de análise de sobrevivência foi apresentado como um modelo alternativo em aplicações de análise de sobrevivência. As vantagens deste são a
flexibilidade em acomodar diferentes formas da função de taxa de falha e a discriminação de modelos.
Por meio de exemplo com dados reais mostrou-se a utilidade desse novo modelo probabilı́stico.
7
Referências Bibliofráficas.
Carrasco, J. M. F., Ortega, E. M. M., Cordeiro, G. M. (2008). A generalized modified Weibull distribution for lifetime modeling. Computational Statistics and Data Analysis, 53, 450-462.
Colosimo, E.A.; Giolo, S.R. Análise de sobrevivência aplicada. São Paulo: Edgard Blucher, 2006.392p.
Mudholkar, G. S., Srivastava, D. K., Kollia, G. D., 1996. A generalization of the Weibull distribution
with application to the analysis of survival data. Journal of American Statistical Association, 91,
1575-1583.
Sen, P.K.; Singer, J.M. Large sample methods in statistics: An introduction with applications.
New York: Chapman and Hall, 1993. 382 p.
Silva, A.N.F. Estudo evolutivo das crianças expostas ao HIV e notificadas pelo núcleo de vigilância
epidemiológica do HCFMRPUSP. 2004. 90 p. Dissertação (Mestrado em Saúde da comunidade).
Universidade de São Paulo, Ribeirão Preto, 2004.
Silva, G. O., Ortega, E. M. M. e Cordeiro, G. M. The Beta Modified Weibull Distribution. Lifetime
Data Analysis, 2010. Aceito para publicação.
6