Obtenir un MTBF ou un indicateur équivalent via une modélisation de Weibull avec censure

Postez ici vos questions, réponses, commentaires ou suggestions - Les sujets seront ultérieurement répartis dans les archives par les modérateurs

Modérateur : Groupe des modérateurs

mostafa benmoussa
Messages : 1
Enregistré le : 29 Nov 2018, 09:38

Obtenir un MTBF ou un indicateur équivalent via une modélisation de Weibull avec censure

Messagepar mostafa benmoussa » 04 Déc 2018, 10:22

Bonjour à tous,

Je me permets de poser une question car je bloque et je n'arrive pas à avancer sur le sujet.
Je vais partir sur un exemple, monter ce que je souhaite obtenir et où ça bloque. Mon code R utilise le package Weibulltools pour estimer les paramètres de la Weibull et traçage du graphique.
En gros, j'ai 50 avions qui volent depuis deux ans avec une date de début précise (que j'aurais choisi, par exemple, janvier 2016) et une date de fin (décembre 2018), ils ont accumulé 150 000 heures de vol et ont eu un total de 250 pannes sur une pièces lambda (une pompe par exemple).
Ces information me permettent de calculer un MTBF (mean time between failure) qui est égal à

Code : Tout sélectionner

Total des heures de vol / Total des pannes
, ce qui donne un MTBF de 600 heures. En théorie, ça n'a de sens que si le taux de défaillance des pièces est constant, ce qui est vrai sur une courte période (6-12 mois). Cette méthode est "facile" car elle ne demande pas trop de calcul.
Dans mon cas, je voudrais obtenir quelques chose de similaire en partant d'une modélisation de Weibull.
J'ai mes données, avec des censures (cas non complet), le MTTF (mean time to failure) est calculable via les paramètres de la Weibull mais pas moyen d'approcher le MTBF. Ma question, est ce possible de le calculer où ca n'a pas de sens ? pas possible ? pour des beta différent de 1 (cas non exponencielle)
Merci bien

CF : http://www.reliawiki.org/index.php/The_Weibull_Distribution

Mon code :

Code : Tout sélectionner

list_of_packages <- c("weibulltools", "lubridate", "dplyr")

new_packages <- list_of_packages[!(list_of_packages %in% installed.packages()[,"Package"])]

if(length(new_packages)){
  install.packages(new_packages)
}

lapply(list_of_packages, require, character.only = TRUE)

#### Simulation weibull vs cumulative MTBF ####

## cass Beta > 1
set.seed(1234)
n <- 1000
b <- 5
e <- 5000

x <- rweibull(n, b, e)

## No censor data, mttf and MTBF
mttf_weib <- round(e * gamma(1 + (1 / b)), 0)
mtbf_cum <- sum(x)/n

mttf_weib
mtbf_cum


## With Censor data
y <- rweibull(n, b, e)

df <- data.frame(x = x, y = y) %>%
  mutate(time = pmin(x, y),
         event = as.numeric(x < y),
         ac = as.character("ouistiti"))

data <- df
# Weibulltools package
distribution <- "weibull"

df_meth <- johnson_method(x = data$time, event = data$event, id = data$ac)
df_meth


df_est <- ml_estimation(x = df_meth$characteristic,
                        event = df_meth$status,
                        distribution = distribution,
                        conf_level = 0.90)

df_est <- rank_regression(x = df_meth$characteristic,
                          y = df_meth$prob,
                          event = df_meth$status,
                          distribution = distribution,
                          conf_level = .95)


conf_betabin <- confint_betabinom(x = df_meth$characteristic,
                                  event = df_meth$status,
                                  loc_sc_params = df_est$loc_sc_coefficients,
                                  distribution = distribution,
                                  bounds = "two_sided",
                                  conf_level = 0.95,
                                  direction = "y")

plot_weibull <- plot_prob(x = df_meth$characteristic,
                          y = df_meth$prob,
                          event = df_meth$status,
                          id = df_meth$id,
                          distribution = distribution,
                          title_main = "Weibull Analysis",
                          title_x = "Flight hours",
                          title_y = "Probability of Failure",
                          title_trace = "Removed Items")

plot_reg_weibull <- plot_mod(p_obj = plot_weibull,
                             x = conf_betabin$characteristic,
                             y = conf_betabin$prob,
                             loc_sc_params = df_est$loc_sc_coefficients,
                             distribution = distribution,
                             title_trace = "Estimated Weibull CDF")

plot_conf_beta <- plot_conf(p_obj = plot_reg_weibull,
                            x = list(conf_betabin$characteristic),
                            y = list(conf_betabin$lower_bound,
                                     conf_betabin$upper_bound),
                            direction = "y",
                            distribution = distribution,
                            title_trace = "Confidence Region")
#### Suppression warnings
plot_conf_beta$elementId <- NULL

ls <- list(plot_weibull = plot_conf_beta,
           n_failure = paste("Nombre de failure" , nrow(data)),
           beta = df_est$coefficients[2],
           eta = df_est$coefficients[1],
           MTTF = round(round(as.numeric(df_est$coefficients[1]), 2) * gamma(1 + (1 / round(as.numeric(df_est$coefficients[2]), 2))), 0),
           data = data)
ls

#Proportion of censors data
table(df$event)
# beta
ls$beta
# eta
ls$eta
# MTBF
ls$MTTF

# 2 times greather ??
sum(df$time)/sum(df$event)

Retourner vers « Questions en cours »

Qui est en ligne

Utilisateurs parcourant ce forum : Aucun utilisateur enregistré et 1 invité