Ejercicio 10.8 Peihua Qiu.

Los datos presentados en la Tabla \(10.4\)son los primeros 20 perfiles observados de un proceso de producción obtenidos durante un monitoreo de proceso de fase II. Los puntos de diseño de todos los perfiles son los mismos, con valores 0.1, 0.25, 0.4, 0.55, 0.7, 0.85, y 1.0. Como en el ejemplo \(10.3\), se cree que el modelo de regresión no lineal de 2 parámetros \((10.30)\) es apropiado para describir estos perfiles. Utilice el gráfico de Shewhart de fase II \((10.26)\) y el gráfico EWMA de fase II \((10.27)–(10.28)\) para monitorear los valores estimados \({(b_{j},d_{j}), j = 1,2,...,20}.\) En ambos gráficos, use \(χ^2_{0.995,2} = 10.597\) como límites de control.

Loas datos de la tabla 10.4 son los sigientes:

Tabla de datos Ejercicio 10.8

Para el menejo de dichos datos hacemos lo sigiente:

ruta<-"C:\\Users\\jdort\\OneDrive - ut.edu.co\\Escritorio\\Base de datos\\20perfilesDatos.xlsx"
library(readxl)
Datos<-read_excel(ruta)
#View(Datos)
attach(Datos) #para tener la disponibilidad de las variables de la tabla de datos

y=cbind(y1j, y2j, y3j, y4j, y5j, y6j, y7j) #une los datos a la derecha de los otros
y
##         y1j   y2j   y3j   y4j   y5j   y6j   y7j
##  [1,] 1.768 1.613 1.421 1.379 1.188 1.113 0.942
##  [2,] 1.890 1.517 1.393 1.299 1.186 1.061 1.074
##  [3,] 1.831 1.597 1.390 1.341 1.085 1.312 0.956
##  [4,] 1.895 1.626 1.506 1.209 1.133 1.009 1.023
##  [5,] 1.702 1.625 1.419 1.466 1.163 1.070 0.931
##  [6,] 1.796 1.618 1.470 1.397 1.273 1.071 1.140
##  [7,] 1.641 1.662 1.376 1.423 1.140 1.213 1.004
##  [8,] 1.630 1.555 1.255 1.308 1.366 0.854 1.098
##  [9,] 1.678 1.782 1.567 1.206 1.150 1.074 0.962
## [10,] 2.076 1.613 1.357 1.354 1.197 1.074 0.991
## [11,] 1.530 1.305 1.177 1.316 0.801 0.072 0.641
## [12,] 1.512 1.513 1.086 1.241 1.002 0.736 0.557
## [13,] 1.556 1.396 1.315 1.233 1.179 0.860 0.694
## [14,] 1.628 1.322 1.177 1.099 1.251 0.788 0.791
## [15,] 1.367 1.294 1.260 1.288 0.960 0.955 0.604
## [16,] 1.445 1.334 1.256 1.276 0.996 0.888 0.775
## [17,] 1.424 1.269 1.260 1.164 1.109 0.845 0.720
## [18,] 1.647 1.334 1.336 1.093 1.048 0.716 0.698
## [19,] 1.457 1.513 1.246 1.181 0.965 0.786 0.819
## [20,] 1.439 1.547 1.337 1.137 1.052 0.867 0.796

Del ejercicio \((10.3)\) se tiene la siguiente informacíón:

El Modelo \((10.30)\) corresponde a \(y=\frac{d}{1+x^{b}}+\epsilon\) para \(x\in[0,1]\) con \(\epsilon\sim N(0,0.25^{2})\) donde \(a\) y \(b\) son dos parámetros con valores \(1\) y \(2\) respectivamente.

Para el desarrollo del ejercicio, se debe tener en cuenta lo siguiente:

Monitoreo de pérfiles no lineal:

  1. Modelo de un perfíl no lineal.

  1. \(y_{ij}=f(x_{i},\theta)+\varepsilon_{ij}\) con \(i=1,2,....,n, ; j=1,2,3,....\) Donde \([{(x_{i},y_{i}), i=1,2,...,n}]\) son \(n\) observaciones de \((x,y)\) , \(f\) es una funcion no lineal conocida, \(\theta\) es un vector de parámetros \(q-dimencional\), y \({\varepsilon_{ij}, i=1,2...,n; j=1,2,3,...}\) son errores aleatorios \(i.i.d.\) con distribución \(N(0,\sigma^2)\)

  2. un modelo podria ser:

    \(f(x_{i}, \theta)=y_{ij}=a+\frac{d}{1+\frac{x_i}{c}}+\varepsilon _{ij}\), \(\theta=(a,c,d)^´\)

  1. Carta \(T^2\)

  1. \(T_{j}^2=\frac{1}{\theta^2}(\hat{\theta_{j}}-\theta_{o})^´(G_{j}^´G_{j})(\hat{\theta_{j}}-\theta_{o})\)

  2. \(\hat{\theta_{j}}=(\hat{\theta_{01}},\hat{\theta_{02}},\hat{\theta_{03}},......\hat{\theta_{0q}})^´\), \(\theta_{0}=(\theta_{01},\theta_{02},\theta_{03},......,\theta_{0q})^´\)

  3. Para un j_ésimo perfil el valor de \(\theta\) puede estimarse a partir de las observaciones \([{(x_{i},y_{ij}),i = 1,2,...,n}]\) utilizando la estimación de mínimos cuadrados no lineales (Para una descricción mas detallada ver el libro de qiu pág-418). Además, el estimador \(\hat{\theta}_{j}\) se puede calcular por la función nls() en R.

  4. \(G_{j}=\left(\frac{\partial f(x_{1},\theta)}{\partial \theta},\frac{\partial f(x_{2},\theta)}{\partial \theta},...,\frac{\partial f(x_{n},\theta)}{\partial \theta} \right)_{n\times q}^´|_{\theta=\theta_{j}}\)

  5. La carta da señal en el momento \(j\) cuando \(T_{j}^2>\chi_{1-\alpha,q}^2\) Donde \(\alpha \in(0,1)\) es un nivel de significación dado, y ,\(q\) es el cuantil \(( 1−\alpha)\)-ésimo de la distribución.

  1. CARTA MEWMA

  1. \(E_{j}=\lambda(\hat{\theta_{j}}-\theta_{o})+(1+\lambda)E_{j-1}\) para \(j\ge 1\) donde \(E_{0} = 0\), y \(\lambda\in(0,1]\) es un parámetro de ponderación.

  2. \(W_{j}^2=E_{j}^´\sum_{E_{j}}^{-1}E_{j}\)

  3. \(\sum_{E_{j}}\approx\lambda^2\sigma^2(G_j^´G_{j})^{-1}+(1-\lambda)^2\sum_{E_{j-1}}\), con \(\sum_{E_{0}}=0\)

  4. La carta produce señal en el tiempo \(j\) si \(W_{j}^2=LCS=\rho\). Donde \(\rho > 0\) es un límite de control elegido para lograr un valor \(ARL_{0}\) dado y \(\sum_{E_{j}}^{}\) es la matriz de covarianza de \(E_{j}\).

Solucion:

##########
#Primero:#  Ajustamos la función NLS para cada perfil... 
##########



x=c(0.1,0.25,0.4,0.55,0.7,0.85,1) #son los puntos de diseño dados.
n = length(x) # son la cantidad de puntos de diseño. 

d = rep(0,20) #crea un vector de ceros de longitud 20
b = rep(0,20) #crea un vector de ceros de longitud 20
mse = array(data=rep(0,80),dim=c(20,2,2)) # crea una matrix con dichas dimenciones
G = matrix(0,n,2) # crea una matrix de ceros de tamaño nx2

for(i in 1:20){
  yy=y[i,]
  fit = nls(yy ~ d1/(1+x^b1),start=list(d1=1,b1=2))
  d[i] = coef(fit)[1]
  b[i] = coef(fit)[2]
  G = cbind(1/(1+x^b[i]),-d[i]*x^b[i]*log(x)/(1+x^b[i])^2)#vector de derivadas                                            evaluado en los coeficientes de cada parámetro.
  mse[i,,] = 0.25^2*solve(t(G)%*%G)
}

#Recordar que la funcion nls() en R determina las estimaciones no lineales #(ponderadas) de mínimos cuadrados de los parámetros de un modelo no lineal.


##########
#Segundo:# se calcula los valores de la estadística de gráficos de Shewhart $T_j^2$
##########

Tj2 = rep(0,20)

for(i in 1:20){
  
  Tj2[i]=t(c(d[i],b[i])-c(2,1))%*%solve(mse[i,,])%*%
    (c(d[i],b[i])-c(2,1))
  
}

##########
#Tercero:# Se calcula los valores de la estadística de gráficos EWMA $E_j$
##########

lambda=0.2
Ej = matrix(0,20,2)
Sigma = array(data=rep(0,80),dim=c(20,2,2))

Ej[1,]=lambda*t(c(d[1],b[1])-c(2,1))
Sigma[1,,]=lambda^2*mse[1,,]
Wj2 = rep(0,20)
Wj2[1] = t(Ej[1,])%*%solve(Sigma[1,,])%*%Ej[1,]

for(i in 2:20){
  
  Ej[i,]=lambda*t(c(d[i],b[i])-c(2,1))+(1-lambda)*Ej[i-1,]
  Sigma[i,,]=lambda^2*mse[i,,]+(1-lambda^2)*Sigma[i-1,,]
  Wj2[i] = t(Ej[i,])%*%solve(Sigma[i,,])%*%Ej[i,]
  
}

#########
#Cuarto:# Se hacen los graficos de las cartas
#########


par(mfrow=c(2,2),mar=c(5,4,1,1)) #se configura el panel del gráfico
j=seq(1,20)


####
#a)# Se grafican los valores estimasdos del parametro b 
####

par(mfrow=c(2,2),mar=c(5,4,1,1))
j=seq(1,20)

plot(j,b,type="p",pch=16,xlab="j",xlim=c(0,21),ylab=expression(hat(b)[j]),
     ylim=c(0,5.5),mgp=c(2,1,0),cex=0.9,xaxt="n",yaxt="n")
axis(1,at=c(0,5,10,15,20),labels=c("0","5","10","15","20"),cex=.7)
axis(2,at=c(0,1,2,3,4,5),labels=c("0","1","2","3","4","5"),cex=.7)
title(xlab="(a)",cex=0.9)

####
#b)# Se grafican los valores estimasdos del parametro d
####

plot(j,d,type="p",pch=16,xlab="j",xlim=c(0,21),ylab=expression(hat(d)[j]),
     ylim=c(0,3),mgp=c(2,1,0),cex=0.9,xaxt="n",yaxt="n")
axis(1,at=c(0,5,10,15,20),labels=c("0","5","10","15","20"),cex=.7)
axis(2,at=c(0,1,2,3),labels=c("0","1","2","3"),cex=.7)
title(xlab="(b)",cex=0.9)


####
#c)# Se grafican la carta shewhart
####


h=qchisq(0.995,2)  ### limite de control para la carta Shewhart

plot(j,Tj2,type="o",lty=1,pch=16,xlab="j",ylab=expression(T[j]^2),
     mgp=c(2,1,0),xlim=c(0,21), ylim=c(0,105),cex=0.9,xaxt="n",yaxt="n")
lines(j,rep(h,20),lty=2,cex=0.9)
axis(1,at=c(0,5,10,15,20),labels=c("0","5","10","15","20"),cex=.7)
axis(2,at=c(0,25,50,75,100),labels=c("0","25","50","75","100"),cex=.7)
title(xlab="(c)",cex=0.9)


####
#d)# Se grafican la carta EWMA 
####

h=qchisq(0.995,2)  ### Limite de control para la carta EWMA

plot(j,Wj2,type="o",lty=1,pch=16,xlab="j",ylab=expression(W[j]^2),
     mgp=c(2,1,0),xlim=c(0,21), ylim=c(0,25),cex=0.9,xaxt="n",yaxt="n")
lines(j,rep(h,20),lty=2,cex=0.9)
axis(1,at=c(0,5,10,15,20),labels=c("0","5","10","15","20"),cex=.7)
axis(2,at=c(0,5,10,15,20,25),labels=c("0","5","10","15","20","25"),cex=.7)
title(xlab="(d)",cex=0.9)

Conclusion:

Se observa que el modelo de regresión no lineal para monitorear los perfiles, en la carta Shewhart presenta un aumento en algunos valores(posiblemente por la naturaleza de los datos), todo lo contrario con respecto a la carta EWMA, la cual se mantiene por debajo del limite de control superior.