Universidade Federal do Rio Grande do Norte
Centro de Ciências Exata e da Terra
Departamento de Estatı́stica
Wanderson Laerte de Oliveira Carvalho
VALIDAÇÃO DE DIAGNÓSTICO PARA
MANUTENÇÃO PREVENTIVA DE POÇOS DE
PETRÓLEO VIA MODELO DE TEMPO DE
FALHA ACELERADA
Natal, dezembro de 2013
Wanderson Laerte de Oliveira Carvalho
VALIDAÇÃO DE DIAGNÓSTICO PARA
MANUTENÇÃO PREVENTIVA DE POÇOS DE
PETRÓLEO VIA MODELO DE TEMPO DE
FALHA ACELERADA
Trabalho apresentado ao Departamento de Estatı́stica
da Universidade Federal do Rio Grande do Norte, em
cumprimento com as exigências legais para obtenção do
tı́tulo de Graduação.
Orientadora:
Profa . Dra . Dione Maria Valença
Natal, dezembro de 2013
ii
Dedicatória
Este trabalho é dedicado único e exclusivamente à pessoa que o estiver lendo. Esta
é a materialização do esforço de um aluno preguiçoso do curso de estatı́stica, que
não acreditava em si, mas com a ajuda de professores, amigos e familiares, obteve a
motivação necessária para concluir esse trabalho.
i
Agradecimentos
Agradeço antes de qualquer coisa aos meus pais e irmãos, pois fizeram de mim quem
eu sou. A minha namorada, foi quem mais me deu apoio moral e me motivou a seguir
em frente, mesmo quando eu estava pronto pra desistir de tudo. Agradeço aos meus
companheiros de curso, por terem me acompanhado em toda a minha trajetória acadêmica, dividindo as cargas e distribuindo motivação. Agradeço aos meus professores,
sem eles eu não teria conseguido o conhecimento necessário para chegar até aqui. Agradeço à minha orientadora Dione Maria Valença pelo suporte e dedicação. Por ultimo,
mas não menos importante agradeço a todos que por mim torceram. Obrigado.
”Exige muito de ti e espera
pouco dos outros. Assim, evitarás
muitos aborrecimentos.”
Confúcio
Resumo
Este trabalho apresenta um estudo em análises de sobrevivência relacionado aos
modelos de tempo de falha acelerada e alguns conceitos sobre técnicas de validação de
testes diagnósticos baseado na curva ROC. O ajuste de um modelo de tempo de falha
acelerado é aplicado nos tempos de falha de poços de petróleo e com base neste modelo
é proposto um teste diagnostico que indica os poços que precisam de manutenção
preventiva. O teste se baseia na sua função de sobrevivência estimada para estimar a
probabilidade de o poço apresentar a primeira falha num intervalo, dado que o poço
funciona no instante t. Por fim, este trabalho avalia o teste proposto baseando-se em
sua sensibilidade e especificidade. O teste necessita de um ponto de corte tal que,
somente os poços que apresentarem probabilidade de falha no intervalo maior que esse
ponto de corte, necessitam de manutenção preventiva. Como resultado, concluı́mos
pela área abaixo da curva ROC que o teste proposto apresenta uma boa capacidade de
diagnosticar poços em iminência de falha.
Palavras-chave: Modelo de tempo de falha acelerada. Poço de petróleo. Validação. Diagnóstico. Curva Roc
iv
Abstract
This work presents a study in survival analyzes related to models of accelerated
failure time and some concepts of techniques about validation of diagnostic tests based
on the ROC curve. The model of accelerated failure time is applied in times of oil
wells fail and based on this model we propose a diagnostic test that indicates the
wells that need preventive maintenance. The test is based on their estimated survival
function to estimate the probability of the well presents the first failure in an interval
fixed, given that it works well at time t. Finally, this paper evaluates the proposed
test based on its sensitivity and specificity. The test requires a cutoff such that only
the wells that show the probability of failure greater than this cut-off range, require
preventive maintenance. As result, we conclude by the area under the ROC curve that
the proposed test has a good ability to diagnose wells on the imminence of failure.
Keywords: Accelerated failure time model. Oil well. Validation. Diagnostic.
ROC curve.
v
Sumário
1 Introdução
1
2 Análise de Sobrevivência
2.1 Função Sobrevivência . . . . . . . . . . . . . . . . . .
2.2 Censura . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 Distribuição Weibull . . . . . . . . . . . . . . . . . .
2.4 Distribuição Valor Extremo . . . . . . . . . . . . . .
2.5 Modelos de Tempo de Falha Acelerada . . . . . . . .
2.6 Inferência em Modelos de Tempo de Falha Acelerada
2.6.1 Estimadores de Máxima Verossimilhança . . .
2.6.2 Testes de Hipóteses . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
3
3
4
5
5
5
6
6
7
.
.
.
.
8
8
9
10
12
4 Aplicação e Discussão
4.1 O Modelo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Aplicação e Validação do Teste Diagnóstico . . . . . . . . . . . . . . . .
14
14
15
5 Considerações finais
21
A Comandos no R
A.1 Leitura de Dados . .
A.2 Divisão da amostra .
A.3 Ajuste do Modelo . .
A.4 Função Sobrevivência
22
22
23
23
24
3 Validação Através da Análise da Curva
3.1 Conceitos Básicos . . . . . . . . . . . .
3.2 A Curva ROC . . . . . . . . . . . . . .
3.3 Interpretação da Curva ROC . . . . . .
3.4 Teste Proposto . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
vi
.
.
.
.
.
.
.
.
ROC
. . . .
. . . .
. . . .
. . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
A.5 Teste Diagnóstico e Função Manut . . . . . . . . . . . . . . . . . . . .
A.6 Gráficos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
A.7 Simulação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
vii
24
25
30
Capı́tulo 1
Introdução
Em analises de sobrevivência existe o interesse de modelar o efeito de variáveis
regressoras no tempo até a ocorrência de um determinado evento de interesse. Os
modelos de tempo de falha acelerada são utilizados em diversas aplicações na literatura.
Ver Colosimo e Giolo (2006).
Nesse trabalho modelamos o efeito de algumas covariáveis no tempo até a ocorrência de falha de equipamentos de sub superfı́cie em poços petrolı́feros, através do modelo
de tempo de falha acelerada, que considera os tempos até a falha dos poços com distribuição Weibull. Baseado no ajusto desse modelo estimamos a probabilidade de um
poço apresentar falha num intervalo de tempo preestabelecido e estabelecemos um teste
diagnóstico que discrimina poços com risco elevado de falha, quando a probabilidade
de falha no intervalo for maior que um determinado valor.
Definido um teste diagnóstico, é importante que este consiga identificar verdadeiramente os indivı́duos portadores da caracterı́stica de interesse, no caso, os poços em
iminência de falha. A análise da curva ROC (Reciver Operating Characteristic) pode
validar o desempenho de um teste diagnóstico, a partir de probabilidades de diagnósticos corretos. O principal objetivo desse trabalho é propor e validar um teste diagnóstico
que classifique com resultado positivo os poços que precisam de manutenção preventiva
e com resultado negativo os que não precisam de manutenção preventiva.
Este trabalho está organizado em 5 capı́tulos, o primeiro já apresentado, sintetiza e
contextualiza toda a informação apresentada. O Capı́tulo 2 apresenta o embasamento
teórico em análise de sobrevivência. São apresentadas algumas distribuições de probabilidade usuais. O modelo de regressão paramétrico capaz de ponderar o efeito de
caracterı́sticas de interesse, no tempo até a falha de unidades experimentais e algumas
noções de inferência em modelos de tempo de falha acelerada.
O Capı́tulo 3 apresenta alguns conceitos básicos de validação de testes diagnósticos
1
2
com base na análise da curva ROC e a formalização do teste diagnóstico proposto.
O Capı́tulo 4 apresenta os resultados da aplicação do teste proposto, num conjunto
de dados reais de tempos até a primeira falha de equipamentos de sub superfı́cie de
poços petrolı́feros.
O Capitulo 5 apresenta as considerações finais do trabalho.
Alem disso, os comandos utilizados no software R são apresentados no Apêndice A
e finalmente as referências.
Capı́tulo 2
Análise de Sobrevivência
A análise de sobrevivência é um dos ramos da estatı́stica que mais vem se destacando. O desenvolvimento de técnicas estatı́sticas somado ao aprimoramento da
capacidade computacional são tidos como motivadores para inúmeras aplicações. Em
análise de sobrevivência a variável resposta é o tempo até ocorrência de um determinado evento. Apesar do termo sobrevivência, nem sempre o evento de interesse é
algo negativo. O tempo estudado pode ser o tempo até a morte de um paciente, mas
também pode ser o tempo até o efeito desejado de uma droga no organismo. Neste
capitulo apresentamos alguns dos conceitos fundamentais em análise de sobrevivência
que são utilizados neste trabalho.
2.1
Função Sobrevivência
Considere T uma variável aleatória que representa o tempo até ocorrência de um
evento de interesse, continua, não negativa e com função densidade de probabilidade
denotada por f (t), onde t é um tempo fixo qualquer. Em análise de sobrevivência é
definida a função de sobrevivência denotada por S(t) e representada pela expressão:
Z∞
S(t) = P (T > t) =
f (t)dt
t
A função sobrevivência deve ser interpretada como a probabilidade de não ocorrência do evento de interesse durante t unidades de tempo. Por exemplo, se o evento
de interesse é a falha de uma pilha, então S(500) equivale à probabilidade de a pilha
não falhar antes de 500 unidades de tempo. Alem disso a função sobrevivência segue
algumas propriedades descritas a seguir:
3
2.2 Censura
4
• S(t) é monótona decrescente
• S(t) é contı́nua a esquerda
• S(0) é 1
• limt7→∞ S(t) = 0
2.2
Censura
Em análise de sobrevivência usualmente o estudo termina antes que todos os elementos apresentem o evento de interesse. Por isso, estudos relacionados à análise de
sobrevivência são caracterizados pela presença de observações incompletas ou parciais,
denominadas censuras, que acontecem por uma variedade de razões, como a não ocorrência do evento de interesse até o termino do estudo ou a perda de acompanhamento
de um elemento estudado. Note que, para indivı́duos censurados, o tempo até ocorrência do evento de interesse, deve ser considerado superior ao tempo registrado de
acompanhamento daquele elemento. Colosimo e Giolo (2006) resaltam o fato de que,
mesmo censurados, todos os resultados provenientes de um estudo de sobrevivência
devem ser usados na análise estatı́stica, e citam duas razões para tal procedimento: (i)
mesmo incompletas, as observações censuradas fornecem informações sobre o tempo de
vida de pacientes; (ii) a omissão das censuras no cálculo das estatı́sticas de interesse
pode acarretar conclusões incorretas.
Existem diferentes tipos de mecanismos de censura, que são regras preestabelecidas
de como será tratada a censura naquele estudo. Os principais tipos são as do tipo
I, do tipo II ou do tipo aleatória. A censura do tipo I é aquela em que o estudo
será terminado após um perı́odo de tempo fixado e todas as observações que ainda
não apresentaram o evento de interesse até esse instante de tempo são censuradas. A
censura do tipo II é aquela em que o estudo será terminado após a ocorrência do evento
de interesse em um numero fixo de observações, depois disso, as observações que ainda
não apresentaram o evento de interesse são consideradas censuras. No estudo aqui
apresentado é utilizado um terceiro mecanismo de censura, a censura do tipo aleatório,
em que a censura é definida quando um individuo, por alguma razão, é retirado do
estudo antes da ocorrência do evento de interesse. Por exemplo, se considerarmos que
o evento de interesses é a morte de um paciente devido ao câncer e o paciente morre
por atropelamento ou outra causa.
2.3 Distribuição Weibull
2.3
5
Distribuição Weibull
Em geral, não se pode assumir que T segue distribuição Normal, pois dados de
sobrevivência apresentam frequentemente distribuição assimétrica positiva. Uma das
distribuições mais comuns para se atribuir a T é a distribuição Weibull.
Uma variável aleatória com distribuição Weibull, com parâmetros α > 0 e γ > 0,
possui funções densidade de probabilidade e de sobrevivência da seguinte forma:
γ
f (t) =
α
γ−1
γ t
t
I(0,∞) (t)
exp −
α
α
γ t
S(t) = exp −
I(0,∞) (t)
α
Uma propriedade interessante dessa distribuição é que a distribuição exponencial pode
ser escrita como uma Weibull, quando o parâmetro γ = 1.
2.4
Distribuição Valor Extremo
Se Y segue distribuição Valor extremo com parâmetros µ e σ, então suas funções
densidade de probabilidade e de sobrevivência são definidas da seguinte forma:
1
f (y) = exp
σ
y−µ
y−µ
−e σ
σ
I(−∞,∞) (y)
n y−µ o
S(y) = exp −e σ I(−∞,∞) (y)
Uma propriedade importante dessa distribuição é que ela pertence à famı́lia de distribuições de posição e escala, ou seja, Y pode ser escrito da forma: Y = µ + σε, em que
ε segue distribuição Valor Extremo Padrão com função densidade:
f0 (ε) = exp {ε − eε }
(2.1)
É possı́vel provar que se T tem distribuição Weibull(α, γ) então log(T ) segue a
distribuição Valor Extremo(µ, σ), com α = eµ e γ = σ1 .
2.5
Modelos de Tempo de Falha Acelerada
Os modelos de tempo de falha acelerada (MTFA) são modelos probabilı́sticos que
procuram representar o efeito das covariáveis no tempo até ocorrência do evento de
2.6 Inferência em Modelos de Tempo de Falha Acelerada
6
Tabela 2.1: Possı́veis Distribuições associadas à Modelos de Tempo de Falha Acelerada
Distribuição de T Distribuição de log T
Exponencial
Valor Extremo
Weibull
Valor Extremo
Log-Normal
Normal
Log-gama
gama
Log-Logistica
Logistica
interesse, ou seja, ele nos indica o quanto aquela caracterı́stica estudada na unidade
amostral, acelera ou retarda o tempo até que ocorra o evento de interesse, o que justifica
o nome do modelo.
Considere T com distribuição tal que log(T ) possui distribuição do tipo posição e
0
escala. Seja x = (1, x1 , x2 , ..., xp ) o vetor de covariáveis que podem afetar o tempo
0
T e seja β = (β0 , β1 , ..., βp ) o vetor de parâmetros associados às covariáveis. Então o
MTFA pode ser escrito da seguinte forma:
log(T ) = x0 β + σε,
ou seja:
log(T ) = β0 + β1 x1 + β2 x2 + ... + βp xp + σε
(2.2)
em que ε tem distribuição padrão que independe de parâmetro desconhecido. A
peculiaridade do modelo nos leva a achar que são poucas as opções de distribuições que
se podem associar a T , mas a verdade é que são várias as possibilidades de escolha, de
tal forma que existem técnicas para obtenção do modelo mais adequada ao conjunto
de dados estudado. A Tabela (2.1) apresenta algumas distribuições que podem ser
atribuı́das à variável T e consequentemente log(T ).
2.6
Inferência em Modelos de Tempo de Falha Acelerada
2.6.1
Estimadores de Máxima Verossimilhança
Ajustar um modelo de tempo de falha acelerada significa estimar os parâmetros
desconhecidos. O método da máxima verossimilhança se torna apropriado em análise de sobrevivência, pois ele incorpora censura, é simples de ser entendido e possui
propriedades ótimas para grandes amostras (Colosimo e Giolo 2006).
2.6 Inferência em Modelos de Tempo de Falha Acelerada
7
Suponha uma amostra não censurada de tempos de falha ti , i ∈ {1, 2, ..., n} de uma
população de interesse, com função de densidade f (ti ; θ), sendo θ vetor de parâmetros.
A função de verossimilhança para θ é definida por:
L(θ) =
n
Y
f (ti , θ)
i=1
.
Segundo o método da máxima verossimilhança, dada uma distribuição, o vetor de
parâmetros que melhor ajusta a distribuição aos dados, é aquele que maximiza a função
de verossimilhança L(θ), que considera que a contribuição de cada ti , em uma amostra
censurada, é a sua função densidade.
Por outro lado, cada observação censurada contribui na verossimilhança com sua
função sobrevivência, o que é apropriado, já que este indivı́duo apresenta tempo de
falha desconhecido, porem, maior que o observado.
Considere o vetor de tempos até a falha ou censura de n unidade observadas e δ o
vetor tal que para cada individuo i δi = 0, se ti é tempo até censura ou δi = 1, se ti é
tempo até falha. Dessa forma fica definida a função de máxima verossimilhança como:
L(θ) =
n
Y
f (ti ; θ)δi S(ti ; θ)1−δi
i=1
2.6.2
Testes de Hipóteses
Muitas vezes há um interesse de testar hipóteses relacionadas aos parâmetros de
um modelo. Por exemplo, no nosso estudo utilizamos testes de hipóteses para testar
se os coeficientes são significativamente diferentes de zero, através do teste de Wald.
Alem disso, a distribuição associada ao modelo ajustado pode ser escolhida através do
teste da razão de verossimilhanças para testar se os parâmetros da Gama Generalizada
assumem valores especı́ficos, uma vez que, dependendo dos valores dos parâmetros a
Gama Generalizada assume a forma de modelos mais simples com: Weibull, Exponencial ou Log-Normal. O Teste da razão de verossimilhanças também é utilizado na
seleção das covariáveis do modelo, em que é testado se a inclusão de cada covariável
altera significativamente o ajuste do modelo. Descrição detalhada da aplicação desses
testes pode ser encontrado em Colosimo e Giolo (2006).
Capı́tulo 3
Validação Através da Análise da
Curva ROC
Pesquisadores frequentemente procuram uma maneira de descrever o quanto um
teste é capaz de diferenciar unidades amostrais quanto a uma determinada caracterı́stica de interesse. Um exemplo simples é o caso dos estudos laboratoriais na medicina,
onde se deseja avaliar o desempenho de um teste diagnóstico para diferenciar doentes
de não doentes.
3.1
Conceitos Básicos
Considere uma amostra de n unidades experimentais providas de uma população
com unidades portadoras e não portadoras de uma determinada caracterı́stica de interesse. Um teste diagnóstico é realizado em cada unidade amostral. Deseja-se que esse
teste apresente resultado positivo quando a unidade é portadora da caracterı́stica de
interesse e negativo, caso contrário. Uma proposta para avaliar esse teste é analisando
as proporções de acerto, que nesse caso é definido de duas formas:
• Acerto positivo: Quando o resultado do teste é positivo e de fato a unidade
observada é portadora da caracterı́stica de interesse.
• Acerto negativo: Quando o resultado do teste é negativo e de fato a unidade
observada não é portadora da caracterı́stica de interesse
Se a avaliação do teste é baseada na proporção de acertos então são definidas duas
formas de avaliação:
8
3.2 A Curva ROC
9
• Sensibilidade (SE ): É a probabilidade do teste fornecer um resultado positivo,
dado que o individuo é realmente portador da caracteristica de interesse. Estimado pela proporção de acertos positivos dentre os portadores da caracterı́stica
de interesse.
• Especificidade (ES ): É a probabilidade do teste fornecer um resultado negativo,
dado que o individuo realmente não porta a caracteristica de interesse. Estimado
pela proporção de acertos negativos dentre os Não portadores da caracterı́stica
de interesse.
Quando os portadores e não portadores da carácteristica de interesse apresentam
uma mesma proporção de acertos positivos, isto é, quando SE é igual a 1 − ES , o teste
é de baixa qualidade, em outras palavras, dado que o teste é positivo, é tão provável
que a unidade testada seja portadora ou não portadora da caracterı́stica de interesse.
O estudo e revisão de medidas usadas para validar o desempenho de testes diagnósticos são dados em Martinez e Louzada-Neto (2000)
3.2
A Curva ROC
De acordo com Matinez et al.(2003), nem sempre o resultado do teste diagnóstico
é dicotômico, muitas vezes esses testes produzem uma variável categórica ordinal ou
contı́nua como resposta. Neste caso, emprega-se uma regra de decisão baseado na
busca de um ponto de corta que discrimina aquela variável resposta como diagnosticada positivamente ou negativamente, tornando a variável resposta dicotômica. Por
exemplo, um indivı́duo com mensurações maiores que o ponto de corte é classificado
como positivo, analogamente, aquele com mensurações menores que o ponto de corte
é classificado como negativo. Desta forma, pode-se calcular SE e ES para diferentes
pontos de corte, dentre os possı́veis valores que o teste produz. Nessa leva de pontos de
corte, um gráfico dos resultantes pares SE e 1 − ES constitui a Curva ROC (Receiver
Operator Characteristic Curve).
Para formalizar esse modelo considere X uma variável aleatória que assume 0 se a
unidade testada não é portadora ou 1 caso contrário e T uma variável aleatória, que
representa o resultado de um teste diagnóstico. Considerando como regra de decisão o
ponto de corte t0 tal que, T > t0 implica na classificação positiva da unidade testada e
T ≤ t0 classifica a unidade testada como negativa. Desta forma, para um especı́fico t0
temos SE = P (T > t0 | X = 1) e ES = P (T ≤ t0 | X = 0). A curva ROC é uma função
contı́nua de SE em 1−ES para uma sequencia de valores t0 contidos no espaço amostral
3.3 Interpretação da Curva ROC
10
de T . Note que, quanto menor o t0 mais facilmente o teste classificará uma unidade
observada como positiva, de tal forma que se todas as unidades observadas forem classificadas como positivas, a proporção de acertos positivos será de 100% (SE = 1), por
outro lado, como todos foram classificado como positivo, a proporção de acertos negativos será de 0% (ES = 0) e estarı́amos classificando como positivo mesmo os negativos
de fato. Muitas vezes o melhor t0 é determinado como aquele que está associado aos
SE e ES simultaneamente maiores, o que nem sempre é adequado. Algumas vezes,
classificar uma unidade como positiva quando ela é na realidade negativa (1 − SE ) é
menos arriscado que classificar como negativa uma unidade que na realidade é positiva
(1 − ES ). Por exemplo, imagine uma vacina que previne uma doença extremamente
letal e que não apresenta efeitos colaterais ou contra indicações, é muito melhor que
o paciente tome a vacina mesmo que nunca venha a ter a doença. Da mesma forma,
existem casos onde 1 − ES pode ser mais arriscado que 1 − SE .
3.3
Interpretação da Curva ROC
Considere os tempos de falha em poços de petróleo t1 , t2 , ..., tn . Dentre estes existem
r poços do tipo A e n−r do tipo B. Considere a variável de interesse X com distribuição
Normal (µA ; σA ) se o poço é do tipo A ou distribuição Normal(µB ; σB ) se é do tipo
B. Não se tem conhecimento se, os poços são da população A ou B, mas são dados as
medidas de X. Considere o teste que indica que um poço é do tipo A quando x for
maior que o ponto de corte x0 . É selecionada uma amostra de poços em que se conhece
o verdadeiro tipo de cada poço, e realiza-se o teste em todos os poços na amostra, para
comparar o resultado do teste com a realidade.
Uma vez definido o teste diagnóstico, estimamos SE = P (X > x0 |A) e ES = P (X ≤
x0 |B) e com auxilio do software R podemos simular as seguintes situações:
• Caso 1: µA próximo de µB . Nesse caso o teste não diferencia bem os poços do
tipo A dos do tpo B, pois a probabilidade de acerto SE é muito próximo da
probabilidade de erro 1 − ES . Essa deficiência é refletida na curva ROC que se
aproxima de uma reta diagonal crescente, indicando que SE é igual a 1 − ES .
• Caso 2: µA um pouco distante de µB . Nesse caso, SE é, em geral, maior que
1 − ES e a curva ROC se aproxima do quadrante superior esquerdo, indicando
que SE é diferente de 1 − ES .
• Caso 3: µA muito distante de µB . Nesse caso, SE fica bem maior que 1 − ES e a
curva ROC apresenta-se ainda mais proximo do quadrante superior esquerdo do
3.3 Interpretação da Curva ROC
11
gráfico.
Assim, quanto mais distantes forem as populações mais fácilmente podemos discriminar seus individuos, e consequentemente melhor os resultados de um teste. Com
relação à curva ROC, temos que quanto maior a capacidade de discriminação do teste,
mais a curva se aproxima do quadrante superior esquerdo. Os gráficos abaixo podem
mostrar de maneira mais clara e simples tudo o que foi exemplificado:
3.4 Teste Proposto
12
Devido a essa relação entre a forma da curva ROC e a qualidade do teste, alguns
autores como Martinez et al. (2003) utilizam a área abaixo da curva ROC como medida
resumo usual do desempenho do teste, resaltando que um teste totalmente incapaz de
discriminar indivı́duos teria uma área sob a curva ROC de 0, 5. Martinez et all (2003)
cita outros autores que consideram a área abaixo da curva ROC como medida não
paramétrica da distancia entre as distribuições dos resultados do teste.
3.4
Teste Proposto
Em nossa aplicação propomos que um teste diagnóstico proposto classifica como
positivo, os poços que precisam de manutenção preventiva e como negativo, os que não
precisam de manutenção preventiva.
Seja T o tempo até a primeira falha de um poço. A regra de classificação é base-
3.4 Teste Proposto
13
ada na probabilidade de o poço apresentar primeira falha em ∆ horas dado que está
funcionando P (t < T < t + ∆|T > t), sendo ∆ a ser fixado pelo setor de manutenção.
O resultado do teste é positivo quando P (t < T < t + ∆|T < t) for maior que uma
determinada probabilidade tolerada de falha (p0 ) e negativo caso contrário.
Considere o modelo (2.4) em que x é o vetor de covariáveis do i-ésimo poço, β o
vetor de coeficientes e σ o parâmetro de escala. Se consideramos com distribuição
Valor Extremo Padrão (ver (2.1)) temos o MTFA Weibull e a função de sobrevivência
de T em função dos parâmetros e do vetor de covariáveis xi dados por:
( σ1 )
t
.
S(t; β, σ, x) = exp − x0 β
e
A probabilidade do poço falhar nas proxima ∆ horas de funcionamento é dada por:
P (t < T ≤ t + ∆|T > t) =
S(t; β, σ, x) + S(t + ∆; β, σ, x)
.
S(t; β, σ, x)
Essa probabilidade está diretamente associada a função de sobrevivência obtida
através do modelo, e portanto a validação do teste também servirá de argumento a
favor da qualidade do modelo.
Ajustado o modelo, tomamos as estimativas de máxima verossimilhança de seus
coeficientes, para estimar a função de sobrevivência em cada poço, da seguinte forma:
( σb1 )
t
bσ
b β,
S(t;
b, x) = exp −
,
ex0 βb
e a partir desta, estimamos as probabilidades de cada poço apresentar a primeira falha
nas proximas ∆ horas de funcionamento:
pb = P (t < T ≤ t + ∆|T > t) =
bσ
bσ
b β,
b + ∆; β,
S(t;
b, x) − S(t
b, x)
.
bσ
b β,
S(t;
b, x)
O resultado do teste é positivo quando essa probabilidade estimada for maior que
um valor de corte pré-estabelecido, indicando que o poço necessita de manutenção
preventiva. Se a probabilidade estimada for menor que o valor de corte, o resultado
do teste será negativo, e o poço é classificado como não necessitado de manutenção
preventiva.
Capı́tulo 4
Aplicação e Discussão
Devido a grande quantidade de poços de petróleo existentes, é de grande importância saber o risco de falha dos poços, com essa informação é desejável indicar os poços
para manutenção antes que falhem. Alem disso, também é de interesse identificar poços
que não necessitam de manutenção, evitando gastos desnecessários e maximizando a
produção de óleo.
4.1
O Modelo
Esse trabalho aborda uma amostra de 603 poços de petróleo, considerando o tempo
até a ocorrência da primeira falha como variável resposta. Algumas caracterı́sticas
dos poços suspeitas de afetar a variável resposta, são incluidas através das covariáveis.
Para analisar a relação entre as covariáveis e o tempo até a falha do poço, é ajustado
um modelo de tempo de falha acelerada Weibull com a seguinte configuração:
log(Ti ) = x0i β + σi
Para i ∈ {1, 2, ..., 603} em que foram consideradas as seguintes covariáveis:
• ProdBase: é a produção base de óleo em metros cúbicos por dia
• ValorBSW: é o percentual de água produzida conjuntamente com o petróleo
• MetEleva: é o método de elevação por bombeio mecânico (BM) ou bombeio por
cavidades progressivas (BCP)
• ValorRGO: é a razão gás-óleo
• IdadePoco: é a idade do poço em anos
14
4.2 Aplicação e Validação do Teste Diagnóstico
15
Tabela 4.1: Resultado do ajuste do modelo de regressão Weibull para uma amostra de
603 poços (dados completos)
Covariável
Valor Erro Padrão Z
P Valor
Intercept
7,148 0,529
13,514 < 0,001
-0,045 0,007
-6,790 < 0,001
ProdBase
0,005 0,002
2,377 0,002
ValorBSW
MetElevaBM
1,772 0,517
3,428 < 0,001
-0,002 0,001
-2,253 0,002
ValorRGO
IdadePoco
0,071 0,008
8,676 < 0,001
0,661 0,329
2,009 0,004
UniOperaOP-CAM
UniOperaOP-ET
2,684 0,792
3,387 < 0,001
1,105 0,475
2,327 0,002
UniOperaOP-RFQ
ProfPoco
0,001 0,001
0,638 0,052
0,002 0,001
3,934 < 0,001
ProfBomb
MEBM:UniOperaOP-CAM 0,834 0,323
2,580 < 0,001
-2,168 0,711
-3,049 < 0,001
MEBM:UniOperaOP-ET
-0,081 0,093
MEBM:UniOperaOP-RFQ -0,038 0,469
MEBM:ProfPoco
-0,002 0,001
-2,968 < 0,001
-1,659 0,009
VBSW:UniOperaOP-CAM -0,005 0,003
VBSW:UniOperaOP-ET
-0,012 0,004
-2,813 < 0,001
-0,576 0,056
VBSW:UniOperaOP-RFQ -0,003 0,005
Log(Escala)
-0,127 0,037
-3,449 < 0,001
• UniOpera: é a unidade operativa do poço. São quatro unidades operacionais:
Alto do Rodrigues (OP-ARG), Canto do Amaro (OP-CAM), Campo de Estreito
(OP-ET) e Riacho da Forquilha (OP-RFQ).
• ProfPoco: e á profundidade do poço em métros
• ProfBomb: é a profundidade em que se encontra instalado a bomba do poço
medido em métros
As estimativas dos coeficientes do ajuste do modelo encontram-se na Tabela 4.1:
Uma descrição mais detalhada do ajuste desse modelo apresentado pode ser encontrada em Dantas et al.(2010).
4.2
Aplicação e Validação do Teste Diagnóstico
Considere a amostra de 603 poços, dividida em duas sub-amostras:
• Amostra de Treinamento (nE = 453): A partir dessa amostra é realizado um novo
ajuste do modelo descrito acima. As estimativas de seus coeficientes encontram-se
4.2 Aplicação e Validação do Teste Diagnóstico
Tabela 4.2: Resultado do ajuste do modelo de
453 poços (amostra de treino)
Covariável
Valor
Intercept
7,45
ProdBase
-0,057
ValorBSW
0,0054
1,47
MetElevaBM
ValorRGO
-0,003
IdadePoco
0,071
0,721
UniOperaOP-CAM
UniOperaOP-ET
2,06
UniOperaOP-RFQ
1,32
ProfPoco
<0,001
0,002
ProfBomb
MEBM:UniOperaOP-CAM 0,863
-1,33
MEBM:UniOperaOP-ET
MEBM:UniOperaOP-RFQ -0,259
MEBM:ProfPoco
-0,002
VBSW:UniOperaOP-CAM -0,006
VBSW:UniOperaOP-ET
-0,014
VBSW:UniOperaOP-RFQ -0,003
-0,121
Log (Escale)
16
regressão Weibull para uma amostra de
Erro Padrão
0,688
0,008
0,003
0,677
0,001
0,009
0,381
0,844
0,633
0,001
0,001
0,378
0,749
0,659
0,001
0,004
0,005
0,005
0,043
Z
10,822
-7,031
2,053
2,167
-2,255
7,302
1,891
2,439
2,086
0.029
3.366
2.281
-1.772
-0.393
-1.741
-1.552
-3.08
-0.54
-2.827
P Valor
< 0,001
< 0,001
0,040
0,030
0,024
< 0,001
0,059
0,015
0,037
0,976
< 0,001
0,023
0,076
0,694
0,082
0,121
0,002
0,589
0,005
4.2 Aplicação e Validação do Teste Diagnóstico
17
na Tabela 4.2:
• Amostra de Teste (nV = 150): Essa amostra de poços será utilizada como suporte no processo de validação, ou seja, os tempos até falha dos poços dessa
amostra serão comparados com os resultados obtidos através do modelo avaliado,
de maneira que, quanto mais os resultados do modelo se aproximarem dos valores
da amostra de teste, melhor será esse modelo. O Gráfico seguinte apresenta a
disperção dos tempos até a primeira falha dos poços da amostra de teste.
Para o processo de validação do modelo são estipulados ∆ ∈ {10.000; 20.000; 30.000}
e p0 ∈ {0, 1; 0, 101; ...; 0, 9}, de maneira que cada ∆ define um teste diferente e quanto
maior o ∆, mais impreciso será o resultado do teste. Para cada p0 é definido as seguintes
formas de acerto, em cada teste:
• Acerto positivo: Quando pb = P (t < T ≤ t + ∆|T > t) >= p0 e o tempo até
falha do i-ésimo poço for menor que ∆, ou seja, quando o resultado do teste for
positivo e aquele poço apresentar falha dentro do intervalo [t; t + ∆].
• Acerto negativo: Quando pb = P (t < T ≤ t + ∆|T > t) < p0 e o tempo até falha
do i-ésimo poço for maior que δ, ou seja, quando o resultado do teste for negativo
e aquele poço apresentar falha fora do intervalo [t; t + ∆].
4.2 Aplicação e Validação do Teste Diagnóstico
18
Os gráficos abaixo mostram as trajetórias das Sensibilidades e Especificidades em
função dos p0 para cada ∆
4.2 Aplicação e Validação do Teste Diagnóstico
19
Esses gráficos podem auxiliar a escolha do melhor valor para p0 , a probabilidade
de falha tolerada, que normalmente é o valor em que as curvas de sensibilidade e
especificidade se cruzam.
Finalmente as curvas ROC para ∆ ∈ {10.000; 20.000; 30.000}, que serviram de
argumento para validação do modelo proposto, são dadas abaixo:
4.2 Aplicação e Validação do Teste Diagnóstico
20
A área abaixo das curvas ROC para ∆ ∈ {10.000; 20.000; 30.000}, calculadas através
da regra do trapézio, são respectivamente: 0.725, 0.750 e 0.883.
Capı́tulo 5
Considerações finais
A curva ROC é um gráfico que cruza a sensibilidade (eixo y) com a inversa da
especificidade (eixo x), para uma sequencia crescente de pontos de corte. Quanto mais
próximo de 1 a área abaixo da curva, melhor a qualidade do teste, por outro lado,
quando mais próximo de 0,5 a área abaixo da curva, pior o teste.
Assim, com base neste modelo ajustado, considerando as caracterı́sticas de cada
poço em funcionamento em um determinado instante, podemos estimar a probabilidade
de falha em um intervalo de ∆ horas, para todos os poços. O teste sugere manutenção
preventiva para aqueles com probabilidade estimada acima de 0,3 se ∆ igual a 10.000
horas e acima de 0.4 para ∆ 20.000 ou 30.000.
A área sob a curva ROC foi de 72%, chegando a 88% quando ∆ = 30.000. Concluı́mos que o teste proposto apresenta uma boa capacidade de diagnosticar poços em
iminência de falha.
21
Apêndice A
Comandos no R
Aqui estão apresentados os comandos utilizados no Software R que resultaram em
toda a extensão deste trabalho
A.1
Leitura de Dados
data = read.table("C:\\Users\\tab2_completa_final.csv" ,sep=";", header = T, dec=",")
# Renomeando as colunas
(Temp = data$tempofalha)
(Cens = data$Indfalha)
(ProdBase = data$CBFP_VL_PROD_BASE)
(ValorBSW = data$CBFP_VL_BSW)
(MetEleva = data$MEEL_CD_METODO)
(ValorRGO = data$CBFP_VL_RGO)
(IdadPoco = data$idade.poco)
(UniOpera = data$UASI_CD_UNID_ADM)
(ProfPoco = data$POCO_VL_PROF_PERF)
(ProfBomb = data$PROF)
(dados = as.data.frame(cbind(Temp,Cens,ProdBase, ValorBSW, MetEleva,
ValorRGO, IdadPoco, UniOpera, ProfPoco, ProfBomb)))
# Classificando os fatores das Covariáveis qualitativas
dados$MetEleva[dados$MetEleva==1] = "BCP"
dados$MetEleva[dados$MetEleva==2] = "BM"
dados$MetEleva<-factor(dados$MetEleva,levels=c("BCP","BM"))
22
A.2 Divisão da amostra
23
dados$UniOpera[dados$UniOpera==1] = "OP-ARG"
dados$UniOpera[dados$UniOpera==2] = "OP-CAM"
dados$UniOpera[dados$UniOpera==3] = "OP-ET"
dados$UniOpera[dados$UniOpera==4] = "OP-RFQ"
dados$UniOpera<-factor(dados$UniOpera,levels=c("OP-ARG","OP-CAM","OP-ET","OP-RFQ"))
summary(dados)
A.2
Divisão da amostra
#Divis~
ao da amostra
id = seq(1:length(dados[,1]))
id.test = sample(id,150)
(dados.test = dados[id.test,])
length(dados.test[,1])
(dados.treino = dados[-id.test,])
length(dados.treino[,1])
A.3
Ajuste do Modelo
#Ajuste do modelodad
require(survival)
ajust1 = survreg(Surv(Temp,Cens) ~ ProdBase + ValorBSW + MetEleva + ValorRGO +
IdadPoco + UniOpera + ProfPoco + ProfBomb + MetEleva*UniOpera +
MetEleva*ProfPoco + ValorBSW*UniOpera,data = dados, dist=’weibull’,x=T)
summary(ajust1) # Original n = 603
ajust2 = survreg(Surv(Temp,Cens) ~ ProdBase + ValorBSW + MetEleva + ValorRGO +
IdadPoco + UniOpera + ProfPoco + ProfBomb + MetEleva*UniOpera +
MetEleva*ProfPoco + ValorBSW*UniOpera,data = dados.treino, dist=’weibull’,x=T)
summary(ajust2)
ajust2.2 = survreg(Surv(Temp,Cens) ~ ProdBase + ValorBSW + MetEleva + ValorRGO +
IdadPoco + UniOpera + ProfBomb + MetEleva*UniOpera +
A.4 Função Sobrevivência
ValorBSW*UniOpera,data = dados.treino, dist=’weibull’,x=T)
summary (ajust2.2)
# Retirada a covariável ProfPoco e a interaç~
ao MetEleva*ProfPoco
ajust3 = survreg(Surv(Temp,Cens) ~ ProdBase + ValorBSW + MetEleva +
ValorRGO + IdadPoco + UniOpera + ProfPoco + ProfBomb + MetEleva*UniOpera +
MetEleva*ProfPoco + ValorBSW*UniOpera,data = dados.test, dist=’weibull’,x=T)
summary(ajust3)
A.4
Função Sobrevivência
#sobreviv^
encia
beta.est = 1/ajust3$scale
alfa.est = exp(ajust3$x %*% ajust2$coef)
sobrev = function(t){
(sobrev = (exp(-(t/alfa.est)^beta.est)))
}
A.5
Teste Diagnóstico e Função Manut
# Teste Diagnóstico
probp = function(t,delt){
(probp = (sobrev(t)-sobrev(t+delt))/sobrev(t))
}
# Fuunç~
ao Manut
manut = function(t0,c0,p0,delt){
cutpoint=NULL
sens=NULL
spec=NULL
invspec=NULL
for(n in seq(from = 0.1, to = 0.9, by = 0.001)){
s=NULL
24
A.6 Gráficos
25
p=p0
t=t0
c=c0
data = as.data.frame(cbind(t,c,p))
data$p[data$p>=n] = 1
# manutenç~
ao
data$p[data$p<n] = 0 # n~
ao manutenç~
ao
data$t[diag(outer((data$t<(min(dados.test$Temp)+delt)), (data$c == 1), "&"))] = 1
data$t[data$t>=(min(dados.test$Temp)+delt)] = 0
data$t[diag(outer((data$t<(min(dados.test$Temp)+delt)), (data$c == 0), "&"))] = NA
data=na.omit(data)
s = data$p+data$t
data = as.data.frame(cbind(t=data$t,p=data$p,s))
cutpoint = c(cutpoint,n)
sens = c(sens,(length(data$s[data$s==2]) / length(data$t[data$t==1])))
spec = c(spec,(length(data$s[data$s==0]) / length(data$t[data$t==0])))
}
return(cbind(cutpoint,sens,spec,invspec=1-spec))
}
A.6
Gráficos
# Curva ROC
# CURVA 1
par(mfrow=c(1,2))
delt = 10000
p1 = c(probp(min(dados.test$Temp),delt))
t1 = dados.test$Temp
c1 = dados.test$Cens
(test1 = manut(t1,c1,p1,delt))
test1 = as.data.frame(test1)
plot(c(0,test1$cutpoint,1), c(1,test1$sens,0),xlab="Ponto de corte",
ylab="Sensibilidade/Specificidade", type="s", pch=20, col=2,
sub="Sensibilidade e Especificidade, Delta = 10000")
A.6 Gráficos
26
points(c(0,test1$cutpoint,1), c(0,test1$spec,1),type="s", col=4, pch=20)
text(0.7,0.9, "Especificidade", cex = 1, col=4)
text(0.7,0.1, "Sensibilidade", cex = 1, col=2)
plot(c(0,1),c(0,1),type="n",axes=1,xlab="1 - Especificidade", ylab="Sensibilidade",
sub="Curva ROC (Delta=10000)")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
#AREA ABAIXO DA CURVA
base=test1[,4]
table(base)
altu=test1[,2]
table(altu)
(B=as.numeric(names(table(base))))
(B=B[-c(7,8,11,14,15,17,19,20,21,27,28,30,31,32,33,35,36,37)])
B=c(0,B,1)
(A=as.numeric(names(table(altu))))
(A=A[c(1,16,18,22,23,24,29,30,33,35,36,37,38,40,45,46,47,48,49,52)])
A=c(0,A,1)
abline(v=B,h=A,col="lightgray")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
soma = 0
for(i in seq(from = 2, to = length(B), by = 1)){
soma = soma + ((B[i]-B[i-1])*A[i-1])
}
soma
#CURVA 2
par(mfrow=c(1,2))
delt = 20000
p1 = c(probp(min(dados.test$Temp),delt))
t1 = dados.test$Temp
c1 = dados.test$Cens
(test1 = manut(t1,c1,p1,delt))
test1 = as.data.frame(test1)
plot(c(0,test1$cutpoint,1), c(1,test1$sens,0),xlab="Ponto de corte",
ylab="Sensibilidade/Specificidade",type="s",pch=20, col=2,
A.6 Gráficos
27
sub="Sensibilidade e Especificidade, Delta = 20000")
points(c(0,test1$cutpoint,1), c(0,test1$spec,1),type= "s", col=4, pch=20)
text(0.7,0.9, "Especificidade", cex = 1, col=4)
text(0.7,0.3, "Sensibilidade", cex = 1, col=2)
plot(c(0,1),c(0,1),type="n",axes=1,xlab="1 - Especificidade", ylab="Sensibilidade",
sub="Curva ROC (Delta=20000)")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
#AREA ABAIXO DA CURVA
base=test1[,4]
table(base)
altu=test1[,2]
table(altu)
(B=as.numeric(names(table(base))))
(B=B[-c(4,10,17,18,19,21,26,27,28,29,31,35)])
B=c(0,B,1)
(A=as.numeric(names(table(altu))))
(A=A[c(6,8,17,25,29,31,34,36,37,38,44,46,48,50,52,55,57,58,59,61,62,63,64)])
A=c(0,A,1)
abline(v=B,h=A,col="lightgray")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
soma = 0
for(i in seq(from = 2, to = length(B), by = 1)){
soma = soma + ((B[i]-B[i-1])*A[i-1])
}
soma
#CURVA3
par(mfrow=c(1,2))
delt = 30000
p1 = c(probp(min(dados.test$Temp),delt))
t1 = dados.test$Temp
c1 = dados.test$Cens
(test1 = manut(t1,c1,p1,delt))
test1 = as.data.frame(test1)
plot(c(0,test1$cutpoint,1), c(1,test1$sens,0),xlab="Ponto de corte",
A.6 Gráficos
28
ylab="Sensibilidade/Specificidade",type="s",pch=20, col=2,
sub="Sensibilidade e Especificidade, Delta = 30000")
points(c(0,test1$cutpoint,1), c(0,test1$spec,1),type="s", col=4, pch=20)
text(0.7,0.9, "Especificidade", cex = 1, col=4)
text(0.7,0.4, "Sensibilidade", cex = 1, col=2)
plot(c(0,1),c(0,1),type="n",axes=1,xlab="1 - Especificidade", ylab="Sensibilidade",
sub="Curva ROC (Delta=30000)")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
#AREA ABAIXO DA CURVA
base=test1[,4]
table(base)
altu=test1[,2]
table(altu)
(B=as.numeric(names(table(base))))
(B=B[-c(5,9,11,13,17)])
B=c(0,B,1)
(A=as.numeric(names(table(altu))))
(A=A[c(24,26,42,43,46,52,55,56,58,59,60,61)])
A=c(0,A,1)
abline(v=B,h=A,col="lightgray")
points(c(1,test1$invspec,0), c(1,test1$sens,0),col=1,type="s",pch=20)
soma = 0
for(i in seq(from = 2, to = length(B), by = 1)){
soma = soma + ((B[i]-B[i-1])*A[i-1])
}
soma
#Gráfico dos tempos
DisTemp = function(t,c){
id = seq(1:length(t))
dados = as.data.frame(cbind(id,t,c))
p1=round(probp(t,10000),2)
p2=round(probp(t,20000),2)
p3=round(probp(t,30000),2)
p=as.data.frame(cbind(p1,p2,p3))
A.6 Gráficos
29
plot(c(0,max(dados$t)+3500),c(0,length(dados$t)),type="n",ylab="Coluna de Produç~
ao",
xlab="Tempo de falha ou censura",axes=1,main="Gráfico de Disperç~
ao dos Tempos")
abline(v=min(dados$t), col="red",lty=2)
#=============================================================== #Os "#" funcionam como interrupto
#===============================================================
abline(v=min(dados$t)+10000, col="red",lty=2)
#abline(v=min(dados$t)+20000, col="red",lty=2)
#abline(v=min(dados$t)+30000, col="red",lty=2)
for(n in seq(from = 1, to = length(t), by = 1)){
if (dados[n,3]==1) {
rect(0,dados[n,1],dados[n,2],dados[n,1],border=1)
points(c(dados[n,2]),c(dados[n,1]),pch=4,col=1)
#===============================================================
# Os "#" funcionam como interruptor deixe sem somente o delta
desejado: 10000,20000 ou 30000 respectivamente
#===============================================================
text(dados[n,2]+3000,dados[n,1],p[n,1], cex =.6)
#text(dados[n,2]+3000,dados[n,1],p[n,2], cex =.6)
#text(dados[n,2]+3000,dados[n,1],p[n,3], cex =.6)
}else{
rect(0,dados[n,1],dados[n,2],dados[n,1],border="grey60",lty=2)
points(c(dados[n,2]),c(dados[n,1]),pch=1,col="grey60")
#===============================================================
# Os "#" funcionam como interruptor deixe sem somente o delta
desejado: 10000,20000 ou 30000 respectivamente
#===============================================================
text(dados[n,2]+3000,dados[n,1],p[n,1], cex =.6,col="grey60")
#text(dados[n,2]+3000,dados[n,1],p[n,2], cex =.6,col="grey60")
#text(dados[n,2]+3000,dados[n,1],p[n,3], cex =.6,col="grey60")
}
}
}
#
#
#
#
id = seq(1:length(dados[,1]))
id.graf = sample(id,25)
(dados.graf = dados[id.graf,])
length(dados.graf[,1])
A.7 Simulação
30
dados.graf = read.table("C:\\Users\\dados.graf.txt", sep=";", header = T)
dados.graf = dados.graf[-1]
t= dados.graf$Temp
c= dados.graf$Cens
DisTemp(t,c)
A.7
Simulação
#SIMULAÇ~
OES
x =
x2=
y =
z =
w =
density(10,
density(11,
density(12,
density(14,
density(10,
bw
bw
bw
bw
bw
=
=
=
=
=
1)
1)
1)
1)
1.05)
Diag = function(liminf,limsup,medianf,mediaf){
SE= NULL
ES= NULL
se= NULL
es= NULL
N = NULL
for(n in seq(from = liminf, to = limsup, by = 0.1)){
se = pnorm(n, mean = mediaf, sd = 1, lower.tail = F, log.p = F)
SE = c(SE,se)
es = pnorm(n, mean = medianf, sd = 1, lower.tail = T, log.p = F)
ES = c(ES,es)
N = c(N,n)
}
return(cbind(N,SE,ES,invES = 1 - ES))
}
#############################################################
dado=as.data.frame(Diag(7,13,10,10))
par(mfrow=c(1,2))
plot(c(7,13),c(0,0.4),type="n",axes=1,main="Curvas de densidade",
xlab="Média",ylab="Densidade")
lines(x,col=4)
A.7 Simulação
lines(w)
abline(v=c(10),lty=2,col="grey")
text(10,0.25, "BCP", cex = 1,col=4)
text(10,0.15, "BM", cex = 1)
text(10,0.05, "to", cex = 1,col="grey")
plot(dado$invES, dado$SE,type="s",main="curva ROC",col=4,
ylab="Sensibilidade",xlab="1 - Especificidade")
####################################################################
dado=as.data.frame(Diag(7,15,10,11))
par(mfrow=c(1,2))
plot(c(7,15),c(0,0.4),type="n",axes=1,main="Curvas de densidade",
xlab="Média",ylab="Densidade")
lines(x,col=4)
lines(x2)
abline(v=c(10,10.5,11),lty=2,col="grey")
text(10,0.2, "BCP", cex = 1,col=4)
text(11,0.2, "BM", cex = 1)
text(10.5,0.1, "to", cex = 1,col="grey")
plot(dado$invES, dado$SE,type="s",main="curva ROC",col=4,
ylab="Sensibilidade",xlab="1 - Especificidade")
####################################################################
dado=as.data.frame(Diag(7,16,10,12))
par(mfrow=c(1,2))
plot(c(7,16),c(0,0.4),type="n",axes=1,main="Curvas de densidade",
xlab="Média",ylab="Densidade")
lines(x,col=4)
lines(y)
abline(v=c(10,11,12),lty=2,col="grey")
text(10,0.2, "BCP", cex = 1,col=4)
text(12,0.2, "BM", cex = 1)
text(11,0.1, "to", cex = 1,col="grey")
plot(dado$invES, dado$SE,type="s",main="curva ROC",col=4,
ylab="Sensibilidade",xlab="1 - Especificidade")
##############################################################
31
A.7 Simulação
dado=as.data.frame(Diag(7,19,10,14))
par(mfrow=c(1,2))
plot(c(7,19),c(0,0.4),type="n",axes=1,main="Curvas de densidade",
xlab="Média",ylab="Densidade")
lines(x,col=4)
lines(z)
abline(v=c(10,12,14),lty=2,col="grey")
text(10,0.2, "BCP", cex = 1,col=4)
text(14,0.2, "BM", cex = 1)
text(12,0.35, "to", cex = 1,col="grey")
plot(dado$invES, dado$SE,type="s",main="curva ROC",col=4,
ylab="Sensibilidade",xlab="1 - Especificidade")
32
Referências
BRAGA, A. C. da S. Curvas ROC: aspectos funcionais e aplicações. 2010. Tese (Doutorado em Métodos Numéricos
e Estatı́sticos) - Universidade do Minho.
COLOSIMO, E.A. e GIOLO, S.R. Análise de Sobrevivência Aplicada. São Paulo : Editora Blucha (Projeto
Fishe,ABE), 2006.
CARNEIRO, H. P. A. Testes de hipóteses em modelos de sobrevivência com fração de cura. 2012. Dissertação
(Mestrado em Matemática Aplicada e Estatı́stica) - Universidade Federal do Rio Grande do Norte.
CARVALHO, J. B. Predição em modelos de tempo de falha acelerado com efeito aleatório para avaliação de riscos
de falha em poços petrolı́feros. 2010. Dissertação (Mestrado em Matemática Aplicada e Estatı́stica) - Universidade
Federal do Rio Grande do Norte.
CARVALHO, J. B.; VALENÇA, D. M.; SINGER, J. M. Prediction of failure probability of oil wells. Submitted to
the Brazilian Journal of Probability and Statistic. 2012.
CASTANHO, M. J. P.; YAMAKAMI, A.; BARROS, L. C.; VENDITE, L. L. Avaliação de um teste em medicina
usando uma Curva ROC Fuzzy. Biomatemática, v. 14, p. 19-28, 2004.
DANTAS, M. A. Modelagem de dados de falhas de equipamentos de sub-superfı́cie em poços de petróleo na bacia potiguar. 2008. Dissertação (Mestrado em Engenharia de Produção) - Universidade Federal do Rio Grande do Norte.
DANTAS, M. A.; VALENÇA, D. M.; FREIRE, M. P. da S.; MEDEIROS, P. G. de; SILVA, D. N. da; ALOISE, D.
J. Modelo de regressão Weibull para estudar dados de falha de equipamentos de sub-superfı́cie em poços petrolı́feros.
Produção (São Paulo. Impresso), v. 20, p. 8, 2010.
ESCOBAR, L. A.; MEEKER, W. Q. A review of accelerated test models. Statistical Science, v. 21, n. 4, p.
552-577, 2006.
FAWCETT, T. An introduction to ROC analysis. Pattern recognition letters, 2006, v. 27, n. 8, p. 861-874.
HULTMANN, R. J. Estudo de confiabilidade de pneus de rodagem. 2011. Dissertação (Programa de Pós-graduação
em Métodos Numéricos em Engenharia) - Universidade Federal do Paraná.
LINNET, K.; BRANDT, E. Assessing diagnostic tests once an optimal cutoff point has been selected. Clinical
chemistry, 1986, v. 32, n. 7, p. 1341-1346.
MARGOTTO, P. R. Curva ROC: como fazer e interpretar no SPSS. Escola Superior de Ciências da Saúde, Distrito
Federal, 2010. Disponı́vel em: <http://www.paulomargotto.com.br/documentos/CurvaR OCS P SS.pdf > .Acessoem :
24abr.2013.
MARTINEZ, E. Z.; LOUZADA-NETO, F.; PEREIRA, B. B. A curva ROC para testes diagnósticos. Cadernos
Saúde Coletiva, v. 11, n. 1, p. 7-31, 2003.
PRATI, R. C.; BATISTA, G. E. A. P. A.; MONARD, M. C. Curvas ROC para a avaliação de classificadores.
Revista IEEE América Latina, v. 6, n. 2, p. 215-222, 2008.
33
A.7 Simulação
34
VAZ, J. C. L. Regiões de incerteza para a Curva ROC em testes diagnósticos. 2009. Dissertação (Mestrado em
Estatı́stica) - Universidade Federal de São Carlos.