forked from flaviobarros/Confiabilidade
-
Notifications
You must be signed in to change notification settings - Fork 0
/
server.R
382 lines (268 loc) · 12.3 KB
/
server.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
# Lógica do servidor
# Inicia o arquivo de servidor com a função principal
shinyServer(function(input, output) {
################## Lendo o conjunto de dados para uso posterior
Data <- reactive({
## Arquivo a ser lido
inFile <- input$file1
## Se nenhum arquivo estiver carregado não faz nada
if (is.null(inFile))
return(NULL)
## Se o arquivo tiver sido carregado, aqui ele é lido
df.raw <- read.csv(inFile$datapath, header=input$header, sep=input$sep, quote=input$quote)
## Cria objeto Surv
## Aqui o padrão é que a primeira coluna tem os tempos
## e a segunda coluna tem o indicador de censura
df.surv <- Surv(df.raw[,1], df.raw[,2])
## Cria o survfit
df.survfit <- survfit(formula = df.surv ~ 1, type="kaplan-meier")
## Cria uma lista para retornar os dados e o Surv
info <- list(df.raw=df.raw, df.surv=df.surv, df.survfit=df.survfit)
return(info)
})
# Opção para tornar paginável
myOptions <- reactive({
list(
page=ifelse(input$pageable==TRUE,'enable','disable'),
pageSize=input$pagesize
)
} )
################## Output para a tabela
output$raw <- renderGvis({
## Por segurança caso o arquivo ainda n tenha sido lido
if (is.null(input$file1)) { return() }
gvisTable(Data()$df.raw, options=myOptions())
})
################## Output para a tabela
output$summary <- renderPrint({
## Segurança para caso nenhum arquivo tenha sido lido
if (is.null(input$file1)) { return() }
## Salvando os dados em um data.frame para ser usado
## localmente nesta função
## Aqui, como a função Data lá em cima para ler, basta
## me referir ao Data() que eu pego o data.frame lido
ys <- Data()$df.survfit
## Apresenta o sumário do modelo de sobrevivência
summary(ys)
})
################## Output do plot 1
output$plot1 <- renderPlot({
## Segurança para caso nenhum arquivo tenha sido lido
if (is.null(input$file1)) { return() }
## Recuperando o survfit
ys <- Data()$df.survfit
## Plot da função de sobrevivência
plot(ys, xlab="Horas", ylab="Probabilidade de sobrevivência")
})
################## Output do plot 2
output$plot2 <- renderPlot({
## Segurança para caso nenhum arquivo tenha sido lido
if (is.null(input$file1)) { return() }
## Recuperando o survfit
dados <- Data()$df.raw
## Pontos de probabilidade igualmente espaçados
p = ppoints(dados$falhas, a=0.3)
## Criando um data.frame com os tempos de falha e o risco acumulado
survresult = data.frame(tempos = dados$falhas, risco = -log(1-p))
## Plot da função de sobrevivência
plot(data = survresult, risco ~ tempos, log="xy", pch=19, col="red", xlab="Horas", ylab="Risco acumulado")
})
################## Output do plot 3
output$plot3 <- renderPlot({
## Segurança para caso nenhum arquivo tenha sido lido
if (is.null(input$file1)) { return() }
## Recebe os dados de sobrevivência
y <- Data()$df.surv
## Estimativa dos parametros da weibull.
##Supondo que os tempos de falha vieram de uma weilbull, então survreg
##encontrar os parametros da weibull que melhor se ajusta aos dados.
##o R calcula estes parametros pelo metodo da maxima verossimilhança
yw = survreg(y ~ 1, dist="weibull")
##Apenas ajustando numero que digitos que se deseja trabalhar.
signif(summary(yw)$loglik, 5)
signif(extractAIC(yw), 5)
## Maximum likelihood estimates:
## For the Weibull model, survreg fits log(T) = log(eta) +
## (1/beta)*log(E), where E has an exponential distribution with mean 1
## eta = Characteristic life (Scale)
## beta = Shape
## Cria lista com estimativas
estimativas <- list()
## Adiciona estimativas a lista
estimativas$escala <- exp(coefficients(yw)[1]) #escala
estimativas$forma <- 1/yw$scale #forma
## Lifetime: expected value and standard deviation.
#\muHAT é o calculo da esperança da minha distribuição
estimativas$média = estimativas$escala * gamma(1 + 1/estimativas$forma)
estimativas$desvio_padrão = estimativas$escala * sqrt(gamma(1+2/estimativas$forma) - (gamma(1+1/estimativas$forma))^2)
## Densidade de probabilidade do modelo ajustado .
curve(dweibull(x, shape=estimativas$forma, scale=estimativas$escala),
from=0, to=estimativas$média+6*estimativas$desvio_padrão, col="blue",
xlab="Horas", ylab="Densidade de probabilidade")
})
################## Output do plot 4
output$plot4 <- renderPlot({
## Criando um gráfico de Cullen and Frey
descdist(Data()$df.raw[,1], boot = 1000)
})
################## Output do plot 5
output$plot5 <- renderPlot({
## Ajustes de distribuição
fw <- fitdist(Data()$df.raw[,1], "weibull")
fg <- fitdist(Data()$df.raw[,1], "exp")
fln <- fitdist(Data()$df.raw[,1], "lnorm")
## Divide a janela em quatro partes
par(mfrow = c(2, 2))
## Cria a legenda
plot.legend <- c("Weibull", "lognormal", "exp")
## Cria os gráficos
denscomp(list(fw, fln, fg), legendtext = plot.legend, xlim = c(0,400), ylim = c(0,0.025), main = 'Histogramas e Densidades Teóricas', xlab = 'dados', ylab = 'Densidade')
qqcomp(list(fw, fln, fg), legendtext = plot.legend, xlab = 'Quantis Teóricos', ylab = 'Quantis Empíricos')
cdfcomp(list(fw, fln, fg), legendtext = plot.legend, main = 'CDFs Empiricas e Teoricas')
ppcomp(list(fw, fln, fg), legendtext = plot.legend, xlab = 'Probabilidade teóricas', ylab = 'Probabilidades empíricas')
})
################## Output dos testes qui-quadrado
output$aderencia <- renderPrint({
## Ajustes de distribuição
fw <- fitdist(Data()$df.raw[,1], "weibull")
fg <- fitdist(Data()$df.raw[,1], "exp")
fln <- fitdist(Data()$df.raw[,1], "lnorm")
## Testes de bondade de ajuste
gofstat(list(fw, fg, fln), fitnames = c('Weibull', 'Exponencial', 'Log-Normal'), discrete = T)
})
################## Output do sumário de Cullen and Frey
output$cullen <- renderPrint({
## Texto de saída do gráfico
descdist(Data()$df.raw[,1], boot = 1000)
})
################## Output das regressão
output$regressao <- renderPrint({
## Recebe os dados de sobrevivência
y <- Data()$df.surv
## Estimativa dos parametros da weibull.
##Supondo que os tempos de falha vieram de uma weilbull, então survreg
##encontrar os parametros da weibull que melhor se ajusta aos dados.
##o R calcula estes parametros pelo metodo da maxima verossimilhança
yw = survreg(y ~ 1, dist="weibull")
## Sumário do ajuste
summary(yw)
})
################## Output das regressão
output$estimativas <- renderPrint({
## Recebe os dados de sobrevivência
y <- Data()$df.surv
## Estimativa dos parametros da weibull.
##Supondo que os tempos de falha vieram de uma weilbull, então survreg
##encontrar os parametros da weibull que melhor se ajusta aos dados.
##o R calcula estes parametros pelo metodo da maxima verossimilhança
yw = survreg(y ~ 1, dist="weibull")
##Apenas ajustando numero que digitos que se deseja trabalhar.
signif(summary(yw)$loglik, 5)
signif(extractAIC(yw), 5)
## Maximum likelihood estimates:
## For the Weibull model, survreg fits log(T) = log(eta) +
## (1/beta)*log(E), where E has an exponential distribution with mean 1
## eta = Characteristic life (Scale)
## beta = Shape
## Cria lista com estimativas
estimativas <- list()
## Adiciona estimativas a lista
estimativas$escala <- exp(coefficients(yw)[1]) #escala
estimativas$forma <- 1/yw$scale #forma
## Lifetime: expected value and standard deviation.
#\muHAT é o calculo da esperança da minha distribuição
estimativas$média = estimativas$escala * gamma(1 + 1/estimativas$forma)
estimativas$desvio_padrão = estimativas$escala * sqrt(gamma(1+2/estimativas$forma) - (gamma(1+1/estimativas$forma))^2)
## Imprime na tela
print(estimativas)
})
################## Output da confiabilidade
output$confiabilidade <- renderPrint({
## Recebe os dados de sobrevivência
y <- Data()$df.surv
## Estimativa dos parametros da weibull.
##Supondo que os tempos de falha vieram de uma weilbull, então survreg
##encontrar os parametros da weibull que melhor se ajusta aos dados.
##o R calcula estes parametros pelo metodo da maxima verossimilhança
yw = survreg(y ~ 1, dist="weibull")
##Apenas ajustando numero que digitos que se deseja trabalhar.
signif(summary(yw)$loglik, 5)
signif(extractAIC(yw), 5)
## Maximum likelihood estimates:
## For the Weibull model, survreg fits log(T) = log(eta) +
## (1/beta)*log(E), where E has an exponential distribution with mean 1
## eta = Characteristic life (Scale)
## beta = Shape
## Cria lista com estimativas
estimativas <- list()
## Adiciona estimativas a lista
escala <- exp(coefficients(yw)[1]) #escala
forma <- 1/yw$scale #forma
## Lifetime: expected value and standard deviation.
#\muHAT é o calculo da esperança da minha distribuição
média = escala * gamma(1 + 1/forma)
desvio_padrão = escala * sqrt(gamma(1+2/forma) - (gamma(1+1/forma))^2)
## Segurança para caso tempo não tenha sido utilizado
if (is.null(input$tempo)) { return() }
## Cálculo da confiabilidade
estimativas$confiabilidade <- round(confiabilidade(forma = forma, escala = escala, tempo = input$tempo),2)
estimativas$confiabilidade <- paste0(estimativas$confiabilidade*100, "%")
## Segurança para caso probabilidade não tenha sido utilizado
if (is.null(input$probabilidade)) { return() }
## Cálculo do quantil
estimativas$quantil <- quantil(escala = escala, forma = forma, probabilidade = 1 - input$probabilidade)
estimativas$quantil <- round(estimativas$quantil, 2)
## Imprime o resultado na tela
print(estimativas)
})
################### Textos nos captions ##################
output$caption1 <- renderText( {
if (is.null(input$file1)) { return() }
"Exploração dos dados"
})
output$caption2 <- renderText( {
if (is.null(input$file1)) { return() }
"Gráfico da função de sobrevivência."
})
output$caption3 <- renderText( {
if (is.null(input$file1)) { return() }
"Sumário do modelo de sobrevivência ajustado"
})
output$caption4 <- renderText( {
if (is.null(input$file1)) { return() }
"Sumário do ajuste para determinação do fator de forma e do fator de escala da Weibull."
})
output$caption5 <- renderText( {
if (is.null(input$file1)) { return() }
"Estimativas de momentos e parâmetros da Weibull."
})
output$caption6 <- renderText( {
if (is.null(input$file1)) { return() }
"Gráfico da função Risco Acumulado."
})
output$caption7 <- renderText( {
if (is.null(input$file1)) { return() }
"Densidade da distribuição Weibull, feito a partir dos parâmetros estimados."
})
output$caption8 <- renderText( {
if (is.null(input$file1)) { return() }
"Calculos da confiabilidade e dos quantis."
})
output$caption9 <- renderText( {
if (is.null(input$file1)) { return() }
"O gráfico de Cullen and Frey fornece a representação gráfica dos momentos de diversas
distribuições de probabilidade e a localização dos dados fornecidos após um processo de reamostragem.
É uma orientação geral para a escolha do melhor modelo probabilístico."
})
output$caption10 <- renderText( {
if (is.null(input$file1)) { return() }
"Os gráficos a seguir apresentam informações com relação a bondade
de ajuste em relação as três distribuições mais utilizadas para modelagem
dos tempos de falha."
})
output$caption11 <- renderText( {
if (is.null(input$file1)) { return() }
"Teste Qui-quadrado para avaliar a bondade de ajuste dos modelos anteriores.
IMPORTANTE: AO TRABALHAR COM DADOS CENSURADOS O PODER DO TESTE DE ADERÊNCIA VAI DIMINUIR."
})
})