Função de intensidade Poisson perturbada pelo número de eventos recorrentes

Função de intensidade Poisson perturbada pelo número de eventos recorrentes

M. Roman^I*; S.L. Caetano^II; F. Louzada^III;J.C. Fogo^I

^IDepartamento de Estatística, DEs-UFSCar, Universidade Federal de São Carlos, Rodovia Washington Luiz, Km 235, 13565-905 São Carlos, SP, Brasil. E-mails: [e-mail]; [e-mail]

IICentro Universitário da Fundação Educacional de Barretos, Avenida Professor Roberto Frade Monte, 389, 14783-226 Barretos, SP, Brasil. E-mail: sabrinalcaetano@gmail

IIIInstituto de Ciências Matemáticas e de Computação, ICMC-USP, Av. Trabalhador Sancarlense, 400, 13566-590 São Carlos, SP, Brasil. E-mail: [e-mail]

RESUMO

Neste trabalho modela-se a função de intensidade de um processo de Poisson considerando o tempo e o total de recorrências, condicionados ao momento anterior, em uma unidade amostral. Adotamos um componente para o processo de Poisson e outro para o número total de eventos ocorridos na unidade. Estudos de simulação e testes de hipóteses empíricos da significância dos parâmetros no modelo foram reali-zados. A significância dos testes de hipótese Wald e da razão de verossimilhança foi aproximadamente 10% para mais de 50 ocorrências. Um conjunto de dados com tempos de recorrência na aquisição de cosméticos foi modelado adequadamente, tendo parâmetros significativos e valores estimados próximos dos valores observados, justificando a utilização do modelo proposto para tempos e números de recorrências em uma unidade amostral.

Palavras-chave: eventos recorrentes, processo de Poisson perturbado, estimação de máxima verossimilhança, teste de Wald, teste de razão de verossimilhança.

ABSTRACT

In this paper we to model the intensity function of Poisson process considering the total number of recurrent events and the total tiem for each sample unit subject to earlier time. Being that one of the component represents the Poisson process and the other the total number of events occurring in the same unit. Simulation studies and empirical hypothesis testing for significance of the parameters were performed. Significance of the Wald and likelihood ratio hypothesis testing were was approximately 10% to more 50 recorrences. A dataset with recorrences time in the acquisition of cosmetics was modeled properly, having significant parameters and estimated and observed values are close, justifying the use of the proposed model for times and numbers of recurrences in a sample unit.

Keywords: Reccurent events, perturbed Poisson process, likelihood estimation, Wald test, likelihood ratio test.

1 INTRODUÇÃO

Uma particularidade da análise de sobrevivência e confiabilidade se refere ao fato de existirem situações nas quais o evento de interesse pode ocorrer repetidas vezes na unidade amostral, que são denominadas dados de eventos recorrentes. Áreas como biomedicina, economia, atuária, criminologia, demografia, engenharia e confiabilidade industrial apresentam dados com essa característica. [5] trazem exemplos específicos nesta área.

Na estrutura de eventos recorrentes podem ser observados o número total de eventos em determinado período de tempo, o tempo exato de recorrência dos eventos, o tempo entre eventos, o custo relacionado a cada recorrência, entre outros [16]. O modelo de sobrevivência pode ser construído considerando um componente aleatório e outro determinístico. O componente aleatório é representado pela distribuição de probabilidade vinculada ao comportamento do tempo de sobrevivência e o componente determinístico é representado pelo relacionamento entre os parâmetros da distribuição e as covariáveis. Os modelos frequentemente utilizados são baseados em processos de Poisson e/ou Renovação, denominados de modelos de intensidade. O processo de Poisson é utilizado quando o interesse é modelar o tempo total de estudo e o processo de Renovação objetiva-se a modelar o tempo intervalar, ver [17].

[11] utilizam processos de Poisson para desenvolver modelos que enfocam o número esperado de eventos ocorridos em determinado intervalo de tempo, [3] sugere modificações nesse modelo. [12] consideram um processo transformado para a função de intensidade acumulada. [4] consideram situações em que, além dos eventos recorrentes, existem eventos terminais que impedem a ocorrência dos demais eventos e, portanto, induzem censuras dependentes.

[19] apresenta um estudo comparativo de modelo de regressão em processos de Poisson. Uma unificação entre processo de Poisson, de Renovação e dos modelos mistos é apresentado por [15], que o denomina de processo de Renovação e Tendência.

[9] optaram por considerar o número total de eventos de recorrência ao invés dos tempos intervalares, ao tratarem o modelo proposto por [12], que traz um modelo com dois processos. Métodos estatísticos para a análise de dados desse tipo não estão completamente disseminados e o assunto ainda é alvo de muitas discussões.

1.1 BACKGROUND

Seja um sistema reparável no qual observamos o tempo t > 0 em que os eventos ocorram nos tempos t₁,t₂,...,t_m, sendo que t₁ < t₂ < ... < t_m < T. Adotamos x_i = t_i-t_i-1 (t₀ = 0 e i = 1,2,3,..., m) para denotar o tempo entre os eventos e N(s,t) para o número de eventos no intervalo [s,t), sendo N(t) = N(0,t) [12]. O modelo de probabilidade para o processo de Poisson é especificado em termos de suas condicionais ou da função de intensidade completa (FIC). Considerando H_t = {N(s):0 < s < t} como o histórico do processo até o instante t. A FIC é dada por

em que H_t = {N(s):0 < s < t} corresponde ao histórico do processo no instante t.

Para este caso, devido a perda de memória, a função intensidade e a função razão são as mesmas, mas isto geralmente não acontece e a função razão não especifica o processo unicamente.

Assumindo um processo de Poisson sendo que ρ(t) é função intensidade (função de risco), N(t) tem distribuição de Poisson com média P(t) = e que o número de eventos recorrentes em cada intervalo de tempo disjunto, são independentes, assim λ(t,H_t) = ρ(t).

Usando a FIC para formular modelos incorporando os tempos de tendência (Poisson) e o número de eventos, geralmente podemos estabelecer efeitos de eventos passados dentro do modelo considerando

em que z(t) = (z₁(t),..., z_p(t))^' é o vetor de funções que podem depender de t e H_t, θ = (θ₁,...,θ_p) é o vetor de parâmetros, ambos desconhecidos.

Considerando (1.1) vários modelos são casos especiais, os quais incluem processos de Poisson com função intensidade dadas por ρ(t) = exp(α + βt), em que z(t) = (1, t)^' e θ = (α, β)^'; ρ(t) = αt^β, em que z(t) = (1, log( t) )^'; e θ = (log (α) , β)^'.

Métodos estatísticos para tais modelos são considerado por vários autores: [1], [2], [6], [8], [10] e [14]. E mesmo não sendo o interesse primordial do trabalho, investigamos os procedimentos de inferência para tais modelos e observamos que as variáveis explicativas podem ser incluídas em termos de z_i(t).

Visando acomodar a influência da quantidade de eventos ocorridos anteriormente condicionado ao momento da ocorrência do evento anterior, propomos um modelo de perturbação no processo de Poisson condicionado.

1.2 Organização do Artigo

Na Seção 2 apresentamos a formulação do modelo utilizando a função intensidade bem como a função de verossimilhança obtida. As estimativas de máxima verossimilhança e a validação do modelo, considerando os testes Wald e Razão de verossimilhança (RV), bem como as estatísticas empíricas associadas a eles, estão apresentadas na Seção 3. Na Seção 4 discutimos o ajuste do modelo proposto para dados artificiais e para dados referentes ao tempo entre compras de produtos de cosméticos. Na Seção 5 alguns comentários finais concluem o artigo.

2 MODELO CONDICIONADO PARA PROCESSO DE POISSON PERTURBADO

Seja ψ um sistema no qual há ocorrência do evento W. Iniciamos o processo no tempo t₀ = 0 e observamos os tempos t_i > 0, i = 1,2,3,...,m, tal que t₁ < t₂ < ... < t_m, nos quais ocorrem o evento de interesse, sendo os tempos t_ii.i.d. Seja N(t_i) o número de eventos ocorridos no intervalo [0, t_i-1], assim, N(t_i) = i-1 é valor relacionado apenas ao tempo presente t_i, sendo constante no tempo anterior, t_i-1.

Consideramos a FIC como apresentada em (1), em que z(t_i) = (1,t_i,i-1)^' e θ = (α, β, γ)^', com α, γ ∈ ℜ e β > 0. A perturbação pelo número de eventos é obtida fazendo g₁(t) = t e g₂(N(t)) = N(t) em

z(t) = (1, g₁(t), g₂[N(t)])',

e a incorporação do efeito dos eventos passados é feito de forma condicional na função intensidade.

A FIC condicionada ao tempo anterior, que expressa o modelo de probabilidade de interesse,é

A expressão (2.1) corresponde ao modelo de intensidade condicional híbrido, em que e^α funciona como a taxa básica, comum para todos os tempos. O parâmetro β é caracterizado como coeficiente da influência do tempo de ocorrência do evento e caracteriza o processo de Poisson. O parâmetro γ indica a influência do número de recorrência no processo e tem significância somente se a função de intensidade depender da quantidade de eventos.

Proposição 2.1.Seja t_i, o tempo de recorrência do evento W, em que t_i-1 < t_i, i = 1,2,...,m, com t₀ = 0, com função de intensidade condicional dada por (2.1). A função densidade de probabilidade condicional do tempo t_i dado t_i-1 é

Proof. Escrevendo f(t_i|t_i-1) = S(t_i|t_i-1)×λ(t_i|t_i-1), com λ(t_i|t_i-1) dada em (2.1), determinamos S(t_i|t_i-1), como segue

portanto

Para λ(u) dada em (2.1), temos

que é uma função de sobrevivência, pois para t_i > 0 temos que

a qual é estritamente negativa, portanto S(t_i|t_i-1) é estritamente decrescente para t_i > 0. E como = 1 e = 0, então F(t_i|t_i-1) = 1 - S(t_i|t_i-1) é uma função densidade de probabilidade acumulada.

Com o produto de (2.3) e (2.1) obtêm-se (2.2) comprovadamente uma f.d.p. □

3 INFERÊNCIA

Assumimos que os tempos t_i, i = 1,2,...,m e m > 0, são i.i.d e realizações da variável aleatória associada ao evento W. Em processos de Poisson condicionados, a função de verossimilhança é dada por L(θ) = λ(t_i|t_i-1)S(t_i|t_i-1), com λ(t_i|t_i-1) dada por (2.1) e S(t_i|t_i-1) dada em (2.3). Sendo θ = (α, β, γ) temos

Para a utilização do método da máxima verossimilhança para a obtenção dos estimadores, utilizamos o logaritmo da função de verossimilhança, a função escore e a matriz de informação observada de Fisher, tudo com respeito a verossimilhança (3.1).

O logaritmo da função verossimilhança representado por ℓ (θ) = log (L (θ)),é

e sua correspondente função escore, obtida da derivada de (3.2) em relação ao vetor paramétrico θ, é dada por

para p = 3, quantidade de parâmetros do vetor θ.

Os elementos da matriz escore, dos quais advêm os estimadores de máxima verossimilhança (EMVs), são

Em princípio se o número total de eventos for grande, então para encontrar as estimativas para θ, a partir de , consideramos que o último tem aproximadamente distribuição Normal com média θ e matriz de variância-covariância I_o()^-1. Escores ou estatísticas de verossimilhança podem ser consideradas e utilizadas de forma usual [7]. Assim, os EMVs não são expressos analiticamente e suas estimativas foram obtidas no software R, pela rotina optim. A obtenção numérica dos EMVs não apresentou problemas, pois a integral envolvida é de fácil resolução, diferentemente do modelo proposto por [12]. A matriz de informação observada I_p×p(θ), que corresponde a derivada de (3.3) em relação ao vetor θ, com r,s = 1,2,..., p, é dada por

Em casos especiais pode-se avaliar a integral em (3.4) analiticamente, mas em geral a integração numérica é necessária [16].

3.1 Consistência dos estimadores de máxima verossimilhança

Fixamos os vetores θ₁ = (2, 0.8, -0.1), θ₂ = (-2, 0.8, -1.5) e θ₃ = (1,0.2, -1.5) de θ = (α, β, γ) e as quantidades de recorrência m = 10, 25, 50,100, 150, 300 e geramos amostras considerando a função inversa da função de probabilidade acumulada condicional. Para avaliar a propriedade de consistência ajustando o modelo aos dados gerados e as estimativas , e , dos respectivos parâmetros, estão apresentadas na Tabela 1. Observamos que os valores estimados se aproximam dos verdadeiros valores dos parâmetros à medida que o tamanho da amostraaumenta. Portanto, segundo [13], os estimadores são consistentes.

3.2 Teste de hipóteses

Utilizamos testes de hipóteses para verificar se existe influência da quantidade de eventos recorrentes no modelo, testando a hipótese H_o: γ = 0. Fixamos o nível de significância de 5% com 999 para θ = θ₁ = (2, 0.8, -0.1) e utilizando os testes Wald e da RV.

Sob a hipótese apresentada, as estatísticas dos testes Wald e da RV são W_γ = e RV = 2l()-2l(0), que assumem assintoticamente a distribuição Qui-Quadrado com 1 grau de liberdade sob H_o, sendo = 3.841 o valor crítico para ambos os testes.

A Tabela 2 dispõem as porcentagens das vezes que o valor da estatística Wald e da estatística RV excederam o valor crítico. Observamos que para 10 recorrências, considerando o teste Wald, 49.6% dos valores são maiores que o valor crítico , contrapondo o esperado nível de significância nominal de 5%. Melhoras ocorrem a partir de 50 recorrências, todavia os valores encontrados estão muito acima do nível de significância nominal de 5% mesmo para amostras de tamanho 300. Situação pior ocorre para o teste RV, no qual para 300 recorrências o percentual de 16.3% dos valores obtidos são maiores que o valor crítico.

Utilizando o método de bootstrap para obter as distribuições assintóticas das estatísticas envolvidas. Consideramos γ = 0 com α = 2 e β = 0.8 para testar a hipótese sobre γ e obtemos as amostras das estimativas dos parâmetros e dos valores das estatísticas empíricas dos testes Wald e da RV. Considerando as estimativas dos parâmetros não fixados pelo teste, α e β, gerou-se nova amostra e desta calculou-se o valor das estatísticas Wald e da RV para 499 repetições deste processo. Dos resultados obtivemos o percentil 95 que foi comparado com a estatística empírica inicialmente calculada. O procedimento de obtenção das 499 estatísticas dos testes, obtenção do percentil e comparação com a estatística empírica foi repetido 399 vezes resultando no percentual de valores da estatística Wald e RV que excederam o valor empírico, cujos resultados apresentamos na Tabela 2, colunas 3 e 5.

Observamos que os valores obtidos, a partir do bootstrap, são próximos do nível de significância, 5% a partir de 50 recorrências. Para a estatística RV, considerando 150 recorrências, essa porcentagem é 6.02% bem próximo do valor esperado. Verificamos que com o aumento do número de recorrência, diminui a porcentagem dos valores que excedem o valor empírico, a partir de um tamanho amostral de 50 observações há uma maior aproximação do valor de significânciade 5% fixado.

4 ANÁLISE DE DADOS

O modelo proposto trata de tempos de recorrência do evento em um mesmo indivíduo ou de um mesmo componente. Nesta seção apresentamos duas análises utilizando o modelo híbrido proposto, inicialmente consideramos dados artificiais e em seguida dados reais referentes aos momentos (tempos) em que um cliente efetua compra de cosméticos. A necessidade de exemplificar o modelo híbrido proposto é de fundamental importância para validar o mesmo em situações em que pode ser utilizado.

4.1 Dados artificiais

Gerou-se 55 recorrências considerando θ = (α, β, γ) = (2, 0.8,-0.1). As EMVs obtidas via maximização direta da função L(θ) são = (2.494,1.081,-0.154) com variância padrão de 0.077, 0.294 e 0.005 para , e respectivamente.

Análise gráfica do ajuste é apresentada na Figura 1(a), a qual expõe o gráfico do número de recorrência versus tempo de recorrência. A linha contínua central apresentada no gráfico é referente ao conjunto de dados obtidos utilizando as estimativas dos parâmetros do modelo ajustado e as demais referem-se aos limitantes de 95% de confiança e os pontos são os valores artificiais utilizados. Assim verificamos que os valores estimados são próximos dos valores artificiais considerados.

O teste de significância dos parâmetros é feita pelo testes de hipóteses Wald, a 5% de significância, cuja estatística é , em que var = (I_3×3 ()^-1) e assume o valor das estimativas , ou , dependendo de qual parâmetro está sendo considerado no teste, e então é comparada com a distribuição Qui-Quadrado e, portanto, o limitante superior da região crítica é = 3.841.

Foram montadas as hipóteses H₀: δ = 0, em que δ = α, β ou γ e as estatísticas calculadas são dadas por = 79.878, = 3.960 e = 4.340, que são maiores que o valor crítico, rejeitando as hipóteses nulas e, portanto, existe risco base e comportamento de tendência nos dados e também há um Processo de Poisson Perturbado, com 95% de confiança.

Para o caso de dados artificiais o modelo teve aplicabilidade e os testes indicam que os parâmetros são significativos ao nível de significância de 5%. Tratamos a seguir um conjunto de dados reais e nele avaliamos o comportamento do modelo.

4.2 Recorrência na compra de cosmético

Na atualidade o setor de cosméticos é um dos segmentos comerciais que mais crescem no Brasil e, nesta seção, tratamos de um exemplo real, que refere-se ao histórico de compras de um cliente de cosméticos. O evento de interesse é o tempo das compras, medido em dias. O conjunto de dados foi obtido junto à uma revendedora e corresponde ao histórico de sete anos de um cliente que teve 57 recorrências. Os dados estão apresentados em dias na Tabela 3.

Para verificar a influência do número de compras anteriores no tempo da compra futura obtemos as estimativas = (2.472,1.329, – 0.167) com variância padrão respectivamente 0.071, 0.406 e 0.006, que, por serem valores baixos validam as estimativas e temos que L () = 67.052.

Com o teste de hipótese com a estatística Wald com nível de significância de 5%, verificamos se o risco base é significativo, se há tendência no processo e se o número de recorrência influencia o processo, utilizando as hipóteses H₀:δ = 0, versus H_a :δ ≠ 0, em que δ = α, β ou γ, dependendo do parâmetro a ser testado.

As estatísticas calculadas são = 86.067, = 4.350 e = 4.648, sendo todas maiores que o valor crítico = 3.841, dando-nos evidências para rejeitar as hipóteses nulas.Podemos afirmar, com 95% de confiança, que existe risco base e comportamento de tendência nos dados, havendo portanto um Processo de Poisson Perturbado.

Com parâmetros significativos, obtivemos o tempo estimado de recorrência, a Figura 1(b) apresenta os valores amostrados e os valores estimados.

Observa-se que o tempo estimado e o tempo observado são próximos, indiferente ao número de compras. Em alguns momentos o valor estimado superestima o valor observado, em outros subestima-o e há situação em que os valores coincidem, o que indica que o modelo proposto é adequado para a análise de eventos que consideram a variável tempo e também o número de eventos ocorridos. Desta forma, estamos satisfeitos com o modelo proposto e a seguir elencamos algumas considerações.

5 COMENTÁRIOS

O modelo proposto é direcionado para análise de tempos de recorrências de um evento na mesma unidade amostral, focando a importância da quantidade de eventos que podem ocorrer na mesma, referindo-se a um Processo de Poisson Perturbado, ou seja, um componente representando o Processo de Poisson e outro o número de eventos ocorridos nesta unidade. De acordo com o estudo de simulação, o modelo proposto é adequado quando trabalhamos com número de recorrências superior a 50. Em relação à consistência dos estimadores, verificamos que quanto maior o número de recorrências mais as estimativas se aproximam do verdadeiro valor do parâmetro. Uma extensão deste estudo é considerar N unidades amostrais e a ocorrência de M_i, i = 1,2,3,...,N recorrências em cada uma delas, tal possibilidade é visionada paratrabalhos futuros.

Recebido em 8 março, 2013

Aceito em 13 dezembro, 2013

• [1] H. Ascher & H. Feingold. "Reparable Systems Reliability: Modelling, Inference, Misconceptions and Their Causes". New York, Marcel Dekker (1984).
• [2] L.J. Bain & M. Engelhardt. Inferences on the parameters and current system reliability for a time truncated Weibull process. Technometrics, 22 (1980), 421-426.
• [3] R. Cook. The design and analysis of randomized trials with recurrent events. Statistics in Medicine, 14 (1995), 2081-2098.
• [4] R. Cook & J. Lawless. Marginal analysis of recurrent events and a terminating event. Statistics in Medicine, 16 (1997), 911-924.
• [5] R.J. Cook & J.F. Lawless. Analysis of repeated events. Statistical Methods in Medical Research 11 (2002), 141-166.
• [6] D.R. Cox & P.A.W. Lewis. "The Statistical Analysis of Series of Events". London: Methuen (1966).
• [7] D.R. Cox & D.V. Hinkley. "Theoretical Statistics". London: Chapman & Hall (1974).
• [8] M.J. Crowder, A.C. Kimber, R.L. Smith & T.J. Sweeting. "Statistical Analysis of Reliability Data". London: Chapman and Hall (1991).
• [9] H. Guo, W. Zhao & A. Mettas. Practical methods for modeling repairable systems with time trends and repair effects. Proceedings of Annual Reliability and Maintainability Symposium, California, (2006), 182-188.
• [10] J.F. Lawless. Regression methods for Poisson process data. Journal of the American StatisticalAssociation, 82 (1987), 807-815.
• [11] J.F. Lawless & C. Nadeau. Some simple robust methods for the analysis of recurrent events. Technometrics, 37 (1995), 158-168.
• [12] J.F. Lawless & K. Thiagarajah. A point-process model incorporating renewals and time trends, with application to repairable systems. Technometrics, 38 (1996), 131-138.
• [13] E.L. Lehmann. "Elements of Large-Sample Theory". New York Springer-Verlag, New York (1999).
• [14] L. Lee & K. Lee. Some results on inference for the Weibull process. Technometrics, 20 (1978), 41-45.
• [15] B.H. Lindqvist, G. Elvebakk & K. Heggland. The trend-renewal process for statistical analysis of repairable systems. Technometrics, 45 (2003), 31-44.
• [16] F. Louzada, J. Mazuchelli & J.A. Achcar. "Introdução à Análise de Sobrevivência e Confiabilidade". III Jornada Regional de Estatística (2002).
• [17] W. Nelson. Confidence Limits for Recurrence Data - Applied to Cost or Number of Product Repair. Technometrics, 37 (1995), 147-157.
• [18] R Development Core Team. R: A Language and Environment for Statistical Computating. R Foundation for Statistical Computing, Vienna, Austria (2011).
• [19] V.L.D. Tomazella. "Modelagem de Dados de Eventos Recorrentes via Processos de Poisson com Termo de Fragilidade", 165p. Tese. (Doutorado em Ciências de Computação e Matemática Computacional) - Instituto de Ciências Matemáticas de São Carlos, Universidade de São Paulo, São Carlos (2003).

*

Autor correspondente: Mari Roman

Datas de Publicação

• Publicação nesta coleção
  07 Mar 2014
• Data do Fascículo
  Dez 2013

Histórico

• Recebido
  08 Mar 2013
• Aceito
  13 Dez 2013