Caso de aplicación

Se estudia la vida útil de la pieza modelo KZY9 que forma parte de un tipo de máquina rotante en una industria, en el que evidentemente la fatiga será en algún grado causa de las fallas. El proveedor envía un set con 250 datos de vida útil, en horas de servicio continuo, de distintas piezas que fueron testeadas en las mismas condiciones de trabajo. Dicha información se encuentra en el archivo de texto adjunto. La planta tiene dos sectores. Se utilizan 50 piezas de este modelo en el sector A y 25 piezas en el sector B. La estrategia de mantenimiento es hacer un recambio de la pieza a las T horas de funcionamiento continuo. Si el recambio se hace antes de la falla tiene un costo C1, mientras que si falla durante la operación el costo es mayor y vale C2. Además, se considera un costo de oportunidad de $3 por cada hora de vida útil no aprovechada para aquellas piezas que sean reemplazadas antes de la falla. Se detallan los valores de costo para los dos sectores de la planta:

SECTOR A SECTOR B
C1 50 50
C2 1000 2000
T 950 750

a)Construya un histograma con los datos de tiempos de falla y calcule el promedio, la varianza muestral y el coeficiente de asimetría. Proponga dos posibles modelos de distribuciones con los que podría describir los tiempos de falla.

Se tiene que \(\overline{X}=\frac{\sum{X_i}}{n}\), \(S^2=\frac{\sum{(X_i - \overline{X})^2}}{n-1}\) y \(\alpha_3 =\frac{E(x-\mu)^3}{(\sqrt{Var(x)})^3}\). Siendo estas las formulas para el cálculo del promedio muestral, la varianza muestral (estimador insesgado) y el coeficiente de asimetria, respectivamente.

A continuación, a partir del ingreso de la muestra, es decir, los datos de la vida útil de 250 piezas modelo KZY9, se obtendra el histograma de la vida útil para este modelo de pieza:

vida_util<-c(908.5, 378.2, 885.2, 1171.9, 625.7, 404.9, 1203.1, 613.5, 955.9, 913.3, 621.6, 735.5, 1248.6, 1054.2, 992.7, 874.4, 722.8, 706.1, 1193.4, 1160.2, 906.6, 1126.8, 605.6, 692.8, 686, 1071.6, 657.3, 1021.4, 750.4, 934.9, 991.9, 1105.2, 898.2, 850.6, 1174.4, 637, 1000.9, 1031.4, 968.6, 301.8, 1124, 1057.1, 901.7, 963.9, 1064.1, 946.2, 538.2, 753.5, 1118.4, 1181.3, 666.4, 1091.1, 736.1, 825.4, 1016.1, 920.3, 841.4, 1181.1, 948.9, 1107.7, 600.9, 1210.4, 1109.8, 1249.3, 1079, 581.8, 558, 1116.1, 948.1, 879, 867.1, 1235.8, 1143.3, 822.6, 1198.1, 626.1, 921.8, 702.1, 1241.7, 794.1, 1140.5, 1021.8, 1079.3, 790.7, 853.3, 1327.7, 1394.7, 1098.1, 1079.6, 820.7, 703.1, 1116.8, 415.2, 857.2, 896.5, 614.8, 1062.3, 316.7, 795.5, 862.5, 681.6, 574, 266.1, 1073, 512.9, 983.3, 1200.8, 760.3, 843.6, 846.6, 897.8, 1119.6, 1039.8, 915.4, 835.9, 751.5, 1165.8, 229, 903.7, 1143.7, 1012.2, 1187.3, 770.6, 726.3, 885.8, 1346.8, 796.6, 599.9, 1045.3, 949.6, 863.5, 1141.2, 870.2, 955.7, 1075.6, 471.9, 1054, 863.6, 865.3, 950.1, 868, 1334, 1347.9, 958.6, 589.3, 573.2, 585.2, 1019.7, 1295.2, 495.3, 1283.3, 724.7, 1215.3, 735.6, 968.5, 820.3, 1052.6, 547.8, 785.1, 1015.5, 1191.2, 1173.1, 708.7, 1085.8, 961.9, 876.8, 941.2, 1013.2, 1226.6, 1145.4, 809.6, 1022.4, 700.5, 1354.1, 1020.1, 1078.3, 711.3, 1197.7, 240.5, 1022.3, 783.7, 816.5, 945.5, 1021.9, 990.3, 853.5, 948.1, 1075.3, 828.4, 755.2, 742.7, 661.4, 1121.5, 998.5, 506.3, 436.6, 1041.1, 572.4, 1050, 1137.5, 1002.9, 1146.9, 888.3, 639, 1333.4, 328.1, 1077.6, 945, 785.1, 1080.5, 1244.3, 789.9, 1278.8, 966.6, 1136.2, 892.6, 1086, 934.4, 857.6, 1115.3, 775.3, 816.7, 1162.8, 981.9, 1175.7, 1088.3, 949.4, 657.7, 621, 1167.5, 901.9, 1101.5, 1076.1, 600.9, 565.1, 1278.9, 1066.5, 701.9, 826.9, 976.3, 875.4, 996.5, 801.5, 910.1, 957.4, 834.1, 935.9, 669.1, 1114.7, 583.9)
hist(vida_util, main = "Histograma de la vida útil de una pieza modelo KZY9", col = "blue", xlab = "Vida útil", ylab = "Frecuencia")

n<-length(vida_util)
prom<-(sum(vida_util)/n)
var<-(sum((vida_util-prom)^2)/(n-1))
desv<-sqrt(var)
library(moments)
asim<-(skewness(vida_util))
tabla<-c(Promedio=prom, Varianza=var, Desvio=desv, Asimetria=asim)
tabla
##      Promedio      Varianza        Desvio     Asimetria 
##   911.0224000 55051.9230304   234.6314621    -0.4905744

Se calcula el promedio muestral, dando un valor igual a 911.022 horas de vida útil, la varianza muestral (insesgada), de valor igual a 55051.92 horas^2 y el coeficiente de asimetria igual a -0.4906. Realizando el Histograma de la vida útil de una pieza KZY9, se puede ver que tiene una asimetria negativa, y corroboramos numericamente esta observacion dado que el coeficiente de asimetria es negativo.

Por lo tanto se proponen dos posibles distribuciones: 1- Weibull: dado que se quiere modelizar el tiempo de vida útil de una pieza que forma parte de un tipo de maquina rotante, entonces el mayor grado de causa de falla en esta pieza sera por fatiga, por lo que se utiliza un modelo de la Weibull de tiempo de aparicion de “fallas por desgaste” (es decir que tendra un parametro w tal que w>1) siendo esta una modelizacion para una distribucion de Weibull de asimetria negativa (es decir con un parámetro w tal que w>3.6) tal como se observa que lo es en el histograma graficado y que se verifica con el valor del coeficiente de asimetria calculado, el cual es negativo. 2- Normal: dado que se toma una muestra n=250 la cual se puede considerar muy grande, pudiendo aproximarce la distribución a una Normal por TCL, ademas podemos observar tanto graficamente en el histograma, como numericamente en valor del coeficiente de asimetria (-0.4906) calculado, que el grado de asimetria negativa es relativamente leve (siendo el valor nulo del coeficiente de asimetria el que indica simetria y no siendo muy lejano al obtenido para este caso).

b)Para los dos modelos de distribución propuestos, estime los parámetros correspondientes usando el método que considere más apropiado en cada caso.

Se realizara una estimación puntual de los respectivos parámetros de las distribuciones propuestas:

1-Weibull: Realizo la estimación puntual por Ajuste de Momentos. Entonces a partir de los datos de la media muestral observada (promedio calculado), \(\hat{\mu}=\overline{x}obs=911.022 horas\), el cual estima el valor de \(\mu\) (media poblacional), y a partrir de la varianza muestral (estimador insesgado) previamente calculada, \(S^2=55051.92 horas^2\), se calcula el desvío muestral (estimador insesgado), dado que \(S=S^2\), el cual estima el valor de \(\sigma\) (desvío poblacional), dando un valor de \(S=234.63 horas\). Luego se obtiene el coeficiente de variación \(\delta=\frac{\sigma}{\mu}=\frac{234.63 horas}{911.022 horas}\) a partir de las estimaciones del desvio y media poblacionales, este da un valor de \(\delta=0.25755\). Entonces voy a la Tabla 4 y (aproximando a un valor de \(\delta=0.25733\)) obtengo el parametro \(\omega=4.4\) (parámetro de forma de la Weibull) y la relación \(\frac{\mu}{\beta}=0.91138\). Habiando obtenido ya la estimación puntual del parametro de forma de la Weibull, a partir de la relación \(\frac{\mu}{\beta}=0.91138\) se obtiene el parametro de escala \(\beta=999.61\), mediante esta estimación puntual por Ajuste de Momentos.

2-Normal: Realizo la estimación puntual, para este caso, los estimadores por Máxima Verosimilitud y por Ajuste de Momentos son los mismos. Entonces habiendo calculado previamente el promedio muestral, se tiene que este es igual a 911.022 horas. Este será el estimador puntual de la media poblacional y a su vez será la estimación del parámetro de posición de la Normal, \(\overline{x}=\hat{\mu}=911.022 horas\). Y luego habiendo calculado tambien la varianza muestral (estimador insesgado), se obtuvo tambien el desvio muestral (estimador insesgado), de un valor de 234.63 horas. Este será el estimador puntual del desvío poblacional y a su vez será la estimación del parámetro de escala de la Normal, \(S=234.63\).

Analizando que tan bien se adaptan estos dos modelos propuestos a la distribución real de la muestra tomada, a continuacion se observan dos “Gráficos Fractiles-Fractiles” (Q-Q Plots) de los valores de la muestra tomada en relación con los valores (fractiles) teóricos tanto de una Weibull (para los parámetros estimados \(\omega=4.4\) y \(\beta=999.61\)) como de una Normal (para los parámetros estimados \(\mu=911.022\) y \(\sigma=234.16\)):

qqplot(rweibull(100000, shape = 4.4, scale = 999.61),vida_util, main = "Gráfico de linealidad Fractiles-Fractiles Weibull", xlab = "Fractiles Teóricos Weibull (simulación de muestra aleatoria de 100000 datos)", ylab = "Fractiles de la muestra", col = "purple")
abline(0,1)

qqplot(rnorm(100000, mean = 911.022, sd = 234.16), vida_util, main = "Gráfico de linealidad Fractiles-Fractiles Normal", xlab = "Fractiles Teóricos Normal (simulación de muestra aleatoria de 100000 datos)", ylab = "Fractiles de la muestra", col = "orange")
abline(0,1)

Entonces cuanto más linealidad se observa en el “gráfico Fractiles-Fractiles”, más se corresponde la distribución de la muestra a la respectiva distribución teórica propuesta. Por lo tanto se observa una levemente mejor linealidad para la distribución Weibull que para la distribución Normal, aunque este es un análisis mas bien cualitativo, por lo que no seria suficiente para decidir aún que distribución sera la mas apropiada para este caso, entonces en el siguiente punto se hará un análisis más cuantitativo para decidir esto.

c)Calcule el logaritmo de la función de verosimilitud (Log-verosimilitud) de la muestra para los dos modelos considerados, y en función del resultado seleccione el mejor modelo.

Sabiendo que la función de verosimilitud viene dada por: \(L(\theta , x)=\prod{f(x_i , \theta)}\) (llamando x a la muestra observada) Entonces el valor de Log-verosimilitud vendra dado por: \(l(\theta , x)= log(L(\theta , x))=\sum{log(f(\theta , x_i))}\)

Para cada caso se calcula el valor de Log-verosimilitud con los estimadores de los respectivos parámetros, previamente calculados:

1-Weibull: \(l(\omega=4.4 , \beta=999.61, x)=\sum{log(f(\omega=4.4, \beta=999.61 , x_i))}= -1714.947\)

log_verosimilitud_weibull<-sum(dweibull(vida_util, shape = 4.4, scale = 999.61, log = TRUE))
tabla_3<-c(LogVerosimilitud=log_verosimilitud_weibull)
tabla_3
## LogVerosimilitud 
##        -1714.947

2-Normal Estandar: \(l(\mu=911.022 , \sigma=234.16, x)=\sum{log(f(\mu=911.022, \sigma=234.63 , x_i))}=-1718.738\)

log_verosimilitud_normal<-sum(dnorm(vida_util, mean = 911.022, sd = 234.63, log = TRUE))
tabla_4<-c(LogVerosimilitud=log_verosimilitud_normal)
tabla_4
## LogVerosimilitud 
##        -1718.739

En función de los resultados obtenidos se llega a la conclusión de que el mejor modelo sera el de la distribución Weibull, dado que la log-verosimilitud de esta es menor que la log-verosimilitud de la distribución Normal Estandar, y por lo tanto la distribución Weibull con los parametros estimados tendra una mayor verosimilitud con la distribución real de la muestra tomada.

d)Para una pieza cualquiera, calcule la media y el desvío estándar del tiempo de servicio.

Habiendo llegado a la conclusión que la distribución que mejor modeliza el tiempo de vida útil de la pieza KZY9 es una distribucion Weibull con parámetros estimados iguales a \(\omega=4.4\) y \(\beta=999.61\), se puede decir que: X:Tiempo de vida útil de un pieza cualquiera KZY9 [en horas].~\(Wei(\omega=4.4 ; \beta=999.61)\)

La media y el desvio de esta variable aleatoria X seran los estimados anteriormente con el promedio muestral y el desvio muestral (estimador insesgado): \(\mu_x=911.022\) y \(\sigma_x=234.63\)

Entonces se define: Y:Tiempo de servicio de una pieza cualquiera KZY9 [en horas].

Se definen los eventos: A:la pieza pertenece al sector A. B:la pieza pertenece al sector B.

Se calcula la probabilidades P(A) y P(B) con la formula de Laplace, sabiendo que hay 50 piezas del sector A y 25 del sector B (75 en total): \(P(A)=\frac{50}{75}=\frac{2}{3}=0.67\) y \(P(B)=\frac{25}{75}=\frac{1}{3}=0.33\)

Se tiene entonces una mezcla de poblaciones (piezas del sector A y piezas del sector B) y por ende se expreza la media de la VA Y con la formula de “esperanza total”: \(E[Y]=E[Y|A]*P(A)+E[Y|B]*P(B)\)

Entonces se definen las VA siguientes: \(Y_A =(Y|A)= X*1[X<T_A]*T_A*1[X>T_A]\) y \(Y_B =(Y|B)= X*1[X<T_B]*T_B*1[X>T_B]\)

T:tiempo de recambio en un determinado sector (no es una VA). De los datos dado: \(T_A =950 horas\) \(T_B =750 horas\)

Se calculan las esperanzas de las dos VA, siendo ambas funciones lineales por tramos de X:

\(E[Y_A]=E[Y|A]=H_w(T_A)+T_A*G_w(T_A)=H_w(950)+950*G_w(950)\)

\(E[Y_B]=E[Y|B]=H_w(T_B)+T_B*G_w(T_B)=H_w(750)+750*G_w(750)\)

Se calculas las respectivas funciones de distribución y esperanzas matemáticas parciales:

\(H_w(950)=\mu*F_\gamma((\frac{950}{\beta})^\omega|1+\frac{1}{\omega} ; 1)=911.022*F_\gamma((\frac{950}{999.61})^(4.4)|1+\frac{1}{4.4} ; 1)=911.022*F_\gamma(0.8|1.23 ; 1)\)

\(H_w(750)=\mu*F_\gamma((\frac{750}{\beta})^\omega|1+\frac{1}{\omega} ; 1)=911.022*F_\gamma((\frac{750}{999.61})^(4.4)|1+\frac{1}{4.4} ; 1)=911.022*F_\gamma(0.28|1.23 ; 1)\)

\(G_w(950)=1-F_w(950|4.4;999.61)\) y \(G_w(750)=1-F_w(750|4.4;999.61)\)

tabla_5<-c(F_gamma_0.8=pgamma(0.8, shape = 1.23, scale = 1), F_gamma_0.28=pgamma(0.28, shape = 1.23, scale = 1), F_weibull_950=pweibull(950, shape = 4.4, scale = 999.61), F_weibull_750=pweibull(750, shape = 4.4, scale = 999.61))
tabla_5
##   F_gamma_0.8  F_gamma_0.28 F_weibull_950 F_weibull_750 
##     0.4472909     0.1603003     0.5503725     0.2461016

Por lo tanto: \(H_w(950)=407.49\) , \(H_w(750)=146.04\) , \(G_w(950)=0.4496\) , \(G_w(750)=0.7539\)

Reemplazando se obtiene la media del tiempo de servicio: \(E[Y]=[H_w(950)+950*G_w(950)]*P(A)+[H_w(750)+750*G_w(750)]*P(B)=[407.49+950*0.4496]*\frac{2}{3}+[146.04+750*0.7539]*\frac{1}{3}\)

=> \(E[Y]=793.56 horas\) Media del tiempo de servicio

Ahora se calcula la varianza de Y por la expresion de la varianza para la mezcla de poblaciones: \(Var[Y]=Var[Y|A]*P(A)+Var[Y|B]*P(B)+P(A)*[E[Y|A]-E[Y]]^2+P(B)*[E[Y|B]-E[Y]]^2\)

De aca se tienen de incognitas las varianzas de (Y|A) y de (Y|B), las cuales se pueden calcular por el Teorema de Steiner: \(Var[Y|A]=Var[Y_A]=E[Y_A ^2]-(E[Y_A])^2\) y \(Var[Y|B]=Var[Y_B]=E[Y_B ^2]-(E[Y_B])^2\)

Definimos las VA siguientes: \(Y_A ^2 =(Y|A)^2= X^2*1[X<T_A]*T_A ^2*1[X>T_A]\) y \(Y_B ^2 =(Y|B)^2= X^2*1[X<T_B]*T_B ^2*1[X>T_B]\)

Entonces se calculan las esperanzas de estas dos VA, siendo ambas funciones lineales por tramos de X:

\(E[Y_A ^2]=H_2w(T_A)+T_A ^2*G_w(T_A)=H_2w(950)+950^2*G_w(950)\)

\(E[Y_B ^2]=H_2w(T_B)+T_B ^2*G_w(T_B)=H_2w(750)+750^2*G_w(750)\)

Habiendo calculado previamente las funciones de distribución necesarias, ahora se calculan las esperanzas matemáticas parciales de segundo orden respectivas:

\(H_2w(950)=(\mu^2+\sigma^2)*F_\gamma((\frac{950}{\beta})^\omega|1+\frac{2}{\omega} ; 1)=(911.022^2+234.63^2)*F_\gamma((\frac{950}{999.61})^(4.4)|1+\frac{2}{4.4} ; 1)=(911.022^2+234.63^2)*F_\gamma(0.8|1.46 ; 1)\)

\(H_2w(750)=(\mu^2+\sigma^2)*F_\gamma((\frac{750}{\beta})^\omega|1+\frac{2}{\omega} ; 1)=(911.022^2+234.63^2)*F_\gamma((\frac{750}{999.61})^(4.4)|1+\frac{2}{4.4} ; 1)=(911.022^2+234.63^2)*F_\gamma(0.28|1.46 ; 1)\)

tabla_6<-c(F_gamma_o.8_alfa1.46=pgamma(0.8, shape = 1.46, scale = 1), F_gamma_o.28_alfa1.46=pgamma(0.28, shape = 1.46, scale = 1))
tabla_6
##  F_gamma_o.8_alfa1.46 F_gamma_o.28_alfa1.46 
##             0.3552932             0.1023961

Por lo tanto: \(H_2w(950)=314438.86\) , \(H_2w(750)=90621.81\)

Reemplazando en Steiner: \(Var[Y|A]=Var[Y_A]=[H_2w(950)+950^2*G_w(950)]-(E[Y_A])^2=[314438.86+950^2*0.4496]-(834.61)^2\)

=> \(Var[Y|A]=Var[Y_A]=23629\)

\(Var[Y|B]=Var[Y_B]=[H_2w(750)+750^2*G_w(750)]-(E[Y_B])^2=[90621.81+750^2*0.7539]-(711.47)^2\)

=> \(Var[Y|B]=Var[Y_B]=8501\)

Ahora se reemplaza todo en la expresión de la varianza para la mezcla de poblaciones:

\(Var[Y]=Var[Y|A]*P(A)+Var[Y|B]*P(B)+P(A)*[E[Y|A]-E[Y]]^2+P(B)*[E[Y|B]-E[Y]]^2\)

\(Var[Y]=23629*0.67+8501*0.33+0.67*[834.61-793.56]^2+0.33*[711.47-793.56]^2\)

=> \(Var[Y]=21989.57horas^2\)

Finalmente se calcula el desvio del tiempo de servicio: => \(\sigma[Y]=148.29 horas\) Desvio del tiempo de servicio

e)En el marco del plan de mantenimiento de esta planta, ¿Cuál es la duración promedio de las piezas que se deben cambiar porque fallan durante la operación?

La duración de las piezas que se debe cambiar porque fallan durante la operación será el tiempo de vida útil de una pieza dado que este tiempo de vida útil sera menor que el tiempo de recambio de la pieza (T). Por lo tanto se define la siguiente VA: (X|X<T): duración de una pieza que se debe cambiar durante la operación.

Se calcula la media de esta VA mediante la expresión de la mezcla de poblaciones (poblaciones de de piezas del sector A y del sector B), es decir haciendo la “esperanza total”: \(E[X|X<T]=E[(X|X<T)|A]*P(A)+E[(X|X<T)|B]*P(B)\)

Se puede expresar del siguiente modo (siendo \(T_A\) y \(T_B\) los tiempos de recambio de los sectores A y B, respectivamente): \(E[X|X<T]=E[X|X<T_A]*P(A)+E[X|X<T_B]*P(B)\)

Ya que estas son esperanzas de VA truncadas, se aplican las siguientes formulas de truncamiento de VA:

\(E[X|X<T_A]=\frac{H_w(T_A)}{F_w(T_A)}\) y \(E[X|X<T_B]=\frac{H_w(T_B)}{F_w(T_B)}\)

Ahora se reemplazan estas expresiónes en el calculo de E[X|X<T], teniendo como datos que \(T_A=950 horas\) y \(T_B=750 horas\): \(E[X|X<T]=\frac{H_w(950)}{F_w(950)}*P(A)+\frac{H_w(750)}{F_w(750)}*P(B)\)

Habiendo calculado estas probabilidades, fuciones de distribución y esperanzas matemáticas parciales previamente (para los parámetros de la Weibull de \(\omega=4.4\) y \(\beta=999.61\)), se reemplazan los valores y se calcula la media: \(E[X|X<T]=\frac{407.49}{0.5504}*0.67+\frac{146.04}{0.2461}*0.33\)

=> \(E[X|X<T]=691.86 horas\) duración promedio de las piezas que se deben cambiar porque fallan durante la operación

f)Para una pieza cualquiera, calcule la media y el desvío estándar del costo de mantenimiento. Analizando el resultado, explique si considera que los valores de T para los dos sectores son apropiados.

Se define la siguiente VA: C:costo de mantenimiento para una pieza cualquiera.

A partir de la fórmula de la media para la mezcla de poblaciones (población del sector A y población del sector B) se tiene: \(E[C]=E[C|A]*P(A)+E[C|B]*P(B)\)

Por lo tanto defino las siguientes VA: \(C_A\)=(C|A):costo de mantenimiento para una pieza del sector A. \(C_A\)=(C|A):costo de mantenimiento para una pieza del sector B.

Considerando los datos dados sobre los costos respectivos para cada sector, en caso de que una pieza falle antes del recambio y en caso de que se le haga el recambio a una pieza antes de que falle, y ademas considerando el costo oportunidad de 3$ por hora de vida útil no aprovechada para las piezas que sean reemplazadas antes de la falla, se tiene que: \(C_A =1000*1[X<T_A]*(50+3*(X-T_A))*1[X>T_A]\) y \(C_B =2000*1[X<T_B]*(50+3*(X-T_B))*1[X>T_B]\)

=> \(C_A =1000*1[X<T_A]*(50+3*X-3*T_A)*1[X>T_A]\) y \(C_B =2000*1[X<T_B]*(50+3*X-3*T_B)*1[X>T_B]\)

Calculo las medias de estas dos VA, siendo ambas fucniones lineales por tramos de X: Para el sector A: \(E[C_A]=1000*F_w(T_A)+50*G_w(T_A)-3*T_A*G_w(T_A)+3*J_w(T_A)\)

=> \(E[C_A]=1000*F_w(T_A)+50*(1-F_w(T_A))-3*T_A*(1-F_w(T_A))+3*(\mu-H_w(T_A))\)

=> \(E[C_A]=1000*F_w(T_A)+50-50*F_w(T_A)-3*T_A+3*T_A*F_w(T_A)+3*\mu-3*H_w(T_A)\)

Con \(T_A=950\) y \(\mu=E[X]=911.022\): \(E[C_A]=1000*F_w(950)+50-50*F_w(950)-3*950+3*950*F_w(950)+3*911.022-3*H_w(950)\)

Con las funciones de distribución y esperanzas parciales previamente calculadas: \(E[C_A]=1000*0.5504+50-50*0.5504-3*950+3*950*0.5504+3*911.022-3*407.49\)

=> \(E[C_A]=802.116\) media del costo de mantenimiento de una pieza del sector A

Para el sector B: \(E[C_B]=2000*F_w(T_B)+50*G_w(T_B)-3*T_B*G_w(T_B)+3*J_w(T_B)\)

=> \(E[C_B]=2000*F_w(T_B)+50*(1-F_w(T_B))-3*T_B*(1-F_w(T_B))+3*(\mu-H_w(T_B))\)

=> \(E[C_B]=2000*F_w(T_B)+50-50*F_w(T_B)-3*T_B+3*T_B*F_w(T_B)+3*\mu-3*H_w(T_B)\)

Con \(T_B=750\) y \(\mu=E[X]=911.022\): \(E[C_B]=2000*F_w(750)+50-50*F_w(750)-3*750+3*750*F_w(750)+3*911.022-3*H_w(750)\)

Con las funciones de distribución y esperanzas parciales previamente calculadas: \(E[C_B]=2000*0.2461+50-50*0.2461-3*750+3*750*0.2461+3*911.022-3*146.04\)

=> \(E[C_B]=1128.566\) media del costo de mantenimiento de una pieza del sector B

Por lo tanto, sabiendo que \(P(A)=0.67\) y \(P(B)=0.33\) y reemplazando los valores obtenidos en la formula para la media de una mezcla de poblaciones, se tiene que: \(E[C]=802.116*0.67+1128.566*0.33\)

=> \(E[C]=909.845\) media del costo de mantenimiento de una pieza cualquiera

Luego para el calculo de la varianza de C, se utiliza la formula para la varianza de una mezcla de poblaciones: \(Var[C]=Var[C|A]*P(A)+Var[C|B]*P(B)+P(A)*(E[C|A]-E[C])^2+P(B)*(E[C|B]-E[C])^2\)

Calculo \(Var[C|A]=Var[C_A]\) y \(Var[C|B]=Var[C_B]\) por Steiner:

\(Var[C_A]=E[C_A^2]-(E[C_A])^2\) y \(Var[C_B]=E[C_B^2]-(E[C_B])^2\)

Entonces se definen las VA \(C_A^2\) y \(C_B^2\) mediante las siguienten expresiones: \(C_A^2 =1000^2*1[X<T_A]*(50+3*(X-T_A))^2*1[X>T_A]\) y \(C_B^2 =2000^2*1[X<T_B]*(50+3*(X-T_B))^2*1[X>T_B]\)

=>\(C_A^2=1000^2*1[X<950]*(2800^2-16800*X+9*X^2)*1[X>950]\) y \(C_B^2=2000^2*1[X<750]*(2200^2-13200*X+9*X^2)*1[X>750]\)

Calculo las medias de estas dos VA, siendo ambas fucniones lineales por tramos de X (siendo \(\mu=E[X]=911.022\) y \(\sigma=\sigma[X]=234.63\)):

Para el sector A: \(E[C_A^2]=1000^2*F_w(950)+2800^2*G_w(950)-16800*J_w(950)+9*J_2w(950)\)

\(E[C_A^2]=1000^2*F_w(950)+2800^2*(1-F_w(950))-16800*(\mu-H_w(950))+9*(\mu^2+\sigma^2-H_2w(950))\)

Reemplazando valores calculados previamente: => \(E[C_A^2]=751087.55\)

Para el sector B: \(E[C_B^2]=2000^2*F_w(750)+2200^2*G_w(750)-13200*J_w(750)+9*J_2w(750)\)

\(E[C_B^2]=2000^2*F_w(750)+2200^2*(1-F_w(750))-13200*(\mu-H_w(750))+9*(\mu^2+\sigma^2-H_2w(750))\)

Reemplazando valores calculados previamente: => \(E[C_B^2]=1685028.202\)

Ahora se calculan las varianzas:

\(Var[C_A]=751087.55-(802.116)^2\) => \(Var[C_A]=107697.47\)

\(Var[C_B]=1685028.202-(1128.566)^2\) => \(Var[C_B]=411366.99\)

Finalmente, reemplazo en la formula de la varianza para la mezcla de poblaciones y calculo el desvio de C: \(Var[C]=Var[C|A]*P(A)+Var[C|B]*P(B)+P(A)*(E[C|A]-E[C])^2+P(B)*(E[C|B]-E[C])^2\)

=> \(Var[C]=107697.47*0.67+411366.99*0.33+0.67*(802.116-909.845)^2+0.33*(1128.566-909.845)^2\)

=> \(Var[C]=231470.95\)

Por lo tanto: \(\sigma[C]=481.11\) desvio del costo de mantenimiento de una pieza cualquiera

Ahora dentro del análisis de los valores de T (tiempo de recambio para cada sector), se calculará los tiempos de recambio óptimos tales que minimizen las respectivas medias de los costos de mantenimiento de una pieza para cada sector. Entonces: Se busca el mínimo de la media de \(C_A\), derivando respecto del tiempo de recmabio para el sector A (siendo este una variable no aleatoria): \(\frac{dE[C_A]}{dT_A}=0=1000*f_w(T_A)-50*f_w(T_A)-3+3*F_w(T_A)+3*T_A*f_w(T_A)-3*T_A*f_w(T_A)\)

=> \(\frac{dE[C_A]}{dT_A}=0=950*f_w(T_A)-3+3*F_w(T_A)\)

=> \(T_A óptimo=906.59horas\) tiempo de recambio óptimo tal que minimiza los costos de mantenimiento del sector A

Del mismo modo se busca el mínimo de la media de \(C_B\), derivando respecto del tiempo de recambio para el secto B (siendo este una variable no aleatoria): \(\frac{dE[C_B]}{dT_B}=0=2000*f_w(T_B)-50*f_w(T_B)-3+3*F_w(T_B)+3*T_B*f_w(T_B)-3*T_B*f_w(T_B)\)

=> \(\frac{dE[C_B]}{dT_B}=0=1950*f_w(T_B)-3+3*F_w(T_B)\)

=> \(T_B óptimo=733.76horas\) tiempo de recambio óptimo tal que minimiza los costos de mantenimiento del sector B

A continuación se presenta una tabla comparativa de los tiempos de recambio óptimos y los tiempo de recambio actuales de los sectores A y B:

SECTOR A SECTOR B
T actual 50 50
T óptimo 1000 2000

Como se puede observar, los tiempos de recambio actuales no son apropiados en comparación a los tiempos de recambio óptimos de la pieza que minimizarían los costos de mantenimiento de cada sector. Por lo tanto, para ambos sectores, se deberían disminuir los tiempos de recambio a los respectivos valores óptimos calculados, \(T_A óptimo=906.59horas\) y \(T_B óptimo=733.76horas\).

g)Se tiene un stock de piezas que ya fueron reemplazadas antes de su vida útil, de las cuales 20 pertenecen al sector A y 10 pertenecen al sector B. ¿Entre qué valores se encontrará, con 90% de probabilidad, la vida útil remanente total de dicho stock?

Para este stock de piezas sucede que el tiempo de recambio fue menor que su vida útil, por ende T<X para estas piezas. Entonces habrá una vida útil remanente, la cual será el tiempo que estas piezas podrían haber seguido en funcionamiento hasta la falla de las mismas. Por lo que esta vida útil remanente vendrá dada por la diferencia entre el tiempo de recambio y la vida útil de la pieza dado que la pieza fue reeemplazada antes de su vida útil. Se definen entonces las siguientes VA:

\(R_i\):vida útil remanente de la i-ésima pieza del stock.

=> \(R_i=(X_i|X_i>T)-T\)

\(R_T\):vida útil remanente total del stock (con un total de 30 piezas en este stock).

=> \(R_T=\sum R_i=\sum[(X_i|X_i>T)-T]\)

Por lo tanto: \((R_i|A)=(X_j|X_j>T_A)-T_A\) y \((R_i|B)=(X_k|X_k>T_B)-T_B\)

Calculo las probabilidades de que una pieza del stock pertenezca al sector A y de que una pieza del stock pertenezca al sector B, mediante Laplace:

\(P(A)=\frac{20}{30}=\frac{2}{3}=0.67\) y \(P(B)=\frac{10}{30}=\frac{1}{3}=0.33\)

Se plantea entonces que la probabilidad de que la vida útil remanente total del stock este en el intervalo [a;b] sea de 0.9 (es decir 90%):

\(P(a<R_T<b)=0.9\) => \(P(a<\sum R_i<b)=0.9\)

Luego si se considera muy grande la cantidad de piezas en este stock (30 piezas en total) y por ende muy grande la cantidad de vidas útiles remanentes sumadas (las cuales son VAs), entonces esta sumatoria, que será igual a la vida útil remanente total, tendrá aproximadamente una distribucion Normal por TCL:

\((R_T=\sum R_i)\)~~\(N(E[R_T];\sigma[R_T])\)

Se calculan entonces la media y la varianza de \(R_T\):

\(E[R_T]=E[\sum R_i]\) => por linealidad => \(E[R_T]=30*E[R_i]\)

\(Var[R_T]=Var[\sum R_i]\) => VAs independiante => \(Var[R_T]=30*Var[R_i]\)

Calculo la media y la varianza de \(R_i\) mediante las respectivas formulas para mezcla de poblaciones: Para la media: \(E[R_i]=E[R_i|A]*P(A)+E[R_i|B]*P(B)\)

=> \(E[R_i]=E[(X_j|X_j>T_A)-T_A]*P(A)+E[(X_k|X_k>T_B)-T_B]*P(B)=(E[X_j|X_j>T_A]-T_A)*P(A)+(E[X_k|X_k>T_B]-T_B)*P(B)\)

Calculo las esperanzas truncada: \(E[X_j|X_j>T_A]=\frac{J_w(T_A)}{G_w(T_A)}=\frac{\mu_x-H_w(950)}{1-F_w(950)}=1119.96\) y \(E[X_k|X_k>T_B]=\frac{J_w(T_B)}{G_w(T_B)}=\frac{\mu_x-H_w(750)}{1-F_w(750)}=1014.7\)

Reemplazando: \(E[R_i]=(1119.96-950)*0.67+(1014.7-750)*0.33\)

=> \(E[R_i]=201.54horas\) media de la vida útil remanente de una pieza cualquiera del stock

=> \(E[R_T]=6046.2horas\) media de la vida útil remanente total del stock

Para la varianza (y el desvio): \(Var[R_i]=Var[R_i|A]*P(A)+Var[R_i|B]*P(B)+P(A)*(E[R_i|A]-E[R_i])^2+P(B)*(E[R_i|B]-E[R_i])^2\)

=> \(Var[R_i]=Var[R_i|A]*P(A)+Var[R_i|B]*P(B)+P(A)*((E[X_j|X_j>T_A]-T_A)-E[R_i])^2+P(B)*((E[X_k|X_k>T_B]-T_B)-E[R_i])^2\)

Calculo las varianzas truncadas por Steiner: \(Var[R_i|A]=Var[(X_j|X_j>T_A)-T_A]=Var[X_j|X_j>T_A]=E[(X_j|X_j>T_A)^2]-(E[X_j|X_j>T_A])^2\)

=>\(Var[R_i|A]=\frac{J_2w(T_A)}{G_w(T_A)}-(E[X_j|X_j>T_A])^2=\frac{\mu_x^2+\sigma_x^2-H_2w(950)}{1-F_w(950)}-1119.96^2\)

=>\(Var[R_i|A]=14758.7\) y \(Var[R_i|B]=Var[(X_k|X_k>T_B)-T_B]=Var[X_k|X_k>T_B]=E[(X_k|X_k>T_B)^2]-(E[X_k|X_k>T_B])^2\)

=>\(Var[R_i|B]=\frac{J_2w(T_B)}{G_w(T_B)}-(E[X_k|X_k>T_B])^2=\frac{\mu_x^2+\sigma_x^2-H_2w(750)}{1-F_w(750)}-1014.7^2\)

=>\(Var[R_i|B]=24092.01\)

Reemplazando: \(Var[R_i]=14758.7*0.67+24092.01*0.33+0.67*((1119.96-950)-201.54)^2+0.33*((1014,7-750)-201.54)^2\)

=>\(Var[R_i]=19864.39horas^2\) varianza de la vida útil remanente de una pieza cualquiera del stock

=>\(Var[R_T]=595931.78horas^2\) varianza de la vida útil remanente total del stock

=>\(\sigma[R_T]=771.97horas\) desvio de la vida útil remanente total del stock

Finalmente usando esta media y desvio calculados, se aplica el TCL y se estandariza la distribució aproximadamente Normal: \(P(\frac{a-\mu[R_T]}{\sigma[R_T]}<Z<\frac{b-\mu[R_T]}{\sigma[R_T]})=0.9\)

Por lo tanto se igualan estos limites inferior y superior a los correspondientes fractiles de la Normal Estandar: \(\frac{a-\mu[R_T]}{\sigma[R_T]}=z0.05=-0.95\) y \(\frac{b-\mu[R_T]}{\sigma[R_T]}=z0.95\)

=> \(a=-z0.95*\sigma[R_T]+\mu[R_T]\) y \(b=z0.95*\sigma[R_T]+\mu[R_T]\)

De este modo,yendo a la tabla 2 se obtiene el fractil \(z0.95=1.6449\) y se reemplaza junto con los valores calculados \(E[R_T]=6046.2horas\) y \(\sigma[R_T]=771.97horas\):

=> \(a=4776.39horas\) y \(b=7316.01horas\)

En conclusión, la vida útil remanente total de este stock estará entre 4776.39 horas y 7316.01 horas, con un 90% de probabilidad.