En un plan de muestreo de aceptación la curva característica de operación muestra la probabilidad de aceptar (\(P_a\)) un lote en función de su nivel de calidad o fracción de disconformidades (\(p\)).
A continuación vamos a ver como representar estas curvas con Julia. También se pueden representar con R, pero en este caso lo más sencillo es utilizar la librería AcceptanceSampling, a no ser que se quieran conocer los detalles de estas representaciones.
La curva característica de operación depende de tres factores:
using Kroki
Diagram(:mermaid, """
flowchart TD
A[Hipergeométrica] --> |"n/N ≤ 0.1 ."| B[p-Binomial]
B --> |"np < 5"| Poisson
B --> |"np ≥ 5"| Normal
A --> |n/N > 0.1 .|C{"¿p?"}
C --> |"p ≤ 0.1 ."| f-Binomial
C --> | p > 0.1 .| Ninguna
""")
Una vez seleccionada la función de distribución a utilizar, solo hay que representar la función de distribución acumulada.
La probabilidad de aceptación es:
\[P_a = \sum_{i=0}^c \frac{\binom{Np}{i} \binom{N-Np}{n-i}}{\binom{N}{n}}\]
donde:
Afortunadamente la librería ‘Distributions.jl’ dispone de la función hipergeométrica y de la función ‘cdf’ para calcular el valor acumulado. La probabilidad de aceptación para un determinado número de unidades disconformes es:
cdf(Hypergeometric(d, N-d, n), c)La representación de la curva característica de operación es sencilla. Consideremos un ejemplo en el que el tamaño del lote es de 100 unidades, el tamaño de la muestra es de 30 y el criterio de aceptación es de 2.
En primer lugar, cargaremos las librerías necesarias y definiremos las variables necesarias:
using Plots, Distributions
N = 100
n = 30
c = 22El paso siguiente será calcular las probabilidades de aceptación para los valores de disconformidades a reprentar. En este caso, realizamos el cálculo para un número de unidades no conformes entre 0 y \(n\):
Pa_h = [cdf(Hypergeometric(d, N-d, n), c) for d in 0:n]31-element Vector{Float64}:
1.0
1.0
1.0
0.9748917748917749
0.9205337617708752
0.8423941179095817
0.7491748936540039
0.64950965921054
0.5504874262796146
0.45720561264903326
0.3728565943582714
0.29905120335385516
0.2362106582040275
⋮
0.03163230398642406
0.02268116633065687
0.016092134445161606
0.011298514978831926
0.00785079636266406
0.005398801057411703
0.0036741919192677004
0.0024744638231675373
0.0016489752345220268
0.001087190671728755
0.000709066446771754
0.000457378125934434Por último solo queda representar la curva característica de operación. Los función de distribución es discreta, por lo solo tenemos los puntos a representar. Añadimos las rectas entre ellos para que la curva quede más clara.
scatter((0:n)./N, Pa_h, legend=false, xlabel="p", ylabel="Pₐ",
title="Función hipergeométrica, N = $N, n = $n, c = $c")
plot!((0:n)./N, Pa_h)En la representación hemos puesto en el eje horizontal \(p\), pero podíamos haber dejado \(d\):
scatter(0:n, Pa_h, legend=false, xlabel="d", ylabel="Pₐ",
title="Función hipergeométrica, n = $n, c = $c")
plot!(0:n, Pa_h)En el caso de utilizar la función de distribución binomial, la función a utilizar es:
\[P_a = \sum_{i=0}^c \binom{n}{i} p^i (1-p)^{n-i}\]
El cálculo de esta probabilidad de aceptación en Julia es muy sencillo. Podemos definir esta función:
Pab(p, n, c) = cdf(Binomial(n, p), c)Pab (generic function with 1 method)Como antes calcularemos las probabilidades de aceptación:
Pa_b = [Pab(p, n, c) for p in (0:n)./N]31-element Vector{Float64}:
1.0
0.9966822906811174
0.9782821654563826
0.9399309054788823
0.8831034382506107
0.8121788131469606
0.7323995229740374
0.6487466279614775
0.5653963645627734
0.4855308229719024
0.41135123955950526
0.34419571712606833
0.2846999592670978
⋮
0.057453385546692645
0.04417898515199697
0.0337085809634231
0.02552420175680366
0.019182250019449253
0.014309325629217927
0.01059587068342088
0.007788748557464795
0.005683539178846932
0.004117070879324837
0.002960509485902006
0.00211317814888654Ya solo queda representar la curva característica de operación de la misma manera que en el caso anterior:
scatter((0:n)./N, Pa_b, legend=false, xlabel="p", ylabel="Pₐ",
title="Función binomial, n = $n, c = $c")
plot!((0:n)./N, Pa_b)En este caso la probabilidad de aceptación se calcula utilizando:
\[P_a = \sum_{i=0}^c \frac{\mu^i \mathrm{e}^{-\mu}}{i!}\]
donde \(\mu\) es el número de defectos medio y su valor es \(\mu = n p\).
Nuevamente definimos una función para calcular la probabilidad de aceptación:
Pap(p, n, c) = cdf(Poisson(p*n), c)Pap (generic function with 1 method)Como en los casos anteriores, calculamos las probabilidades de aceptación a representar y realizamos el gráfico:
Pa_p = [Pap(p, n, c) for p in (0:n)./N]
scatter((0:n)./N, Pa_p, legend=false, xlabel="p", ylabel="Pₐ",
title="Función Poisson, N = $N, n = $n, c = $c")
plot!((0:n)./N, Pa_p)Podemos comparar las tres curvas características de operación:
plot((0:n)./N, [Pa_h, Pa_b, Pa_p], xlabel="p", ylabel="Pₐ",
title="N = $N, n = $n, c = $c",
label=["Hipergeomética" "Binomial" "Poisson"])En el caso de los planes de aceptación para variables utilizaremos la función de distribución normal. A diferencia de las anteriores se trata de una función continua.
Si se conoce la desviación estándar histórica del parámetro de calidad, la probabilidad de aceptación se calcula en este caso como:
\[P_a = \int_Q^\infty \frac{1}{\sqrt{2 \pi}} \mathrm{e}^{-\frac{t^2}{2}}\mathrm{d}t\]
donde \(Q = k+\Phi^{-1}(p)\sqrt(n)\). \(\Phi^{-1}\) es la función quantil o función distribución normal inversa.
Calcularemos la probabilidad de aceptación con la siguente función:
Pan(p, n, k) = ccdf(Normal(), (k+quantile(Normal(),p))*sqrt(n))Pan (generic function with 1 method)Ya estamos en condiciones de representar la curva característica de operación:
k = 1.5
Pa_n = [Pan(p, n, k) for p in 0:0.005:0.2]41-element Vector{Float64}:
1.0
0.9999999980987077
0.9999969958027363
0.9998788367564669
0.9987893528726985
0.9941210214350168
0.9814967563928896
0.9562194005905845
0.9151337992857762
0.857743753688643
0.7862264268854942
0.7046514020040986
0.6179142189650453
⋮
0.007754545121397837
0.005557555519595953
0.00396262028417792
0.002811706079611004
0.001985874005899692
0.001396440504524963
0.0009778385404111605
0.0006819651249039486
0.0004737792691404002
0.00032792188789803315
0.00022615204871624653
0.00015542327098513168plot(0:0.005:0.2, Pa_n, legend=false, xlabel="p", ylabel="Pₐ",
title="σ conocida, n = $n, k = $k")En el caso de de no conocer la desviación estándar del parámetro de calidad, el cálculo de la probabilidad de aceptación es más complejo:
\[P\left(f, \delta, t\right)=\frac{f !}{2^{(f-1) / 2} \Gamma(f / 2) \sqrt{\pi f}} \int_{-\infty}^{t} \mathrm{e}^{-(1 / 2)\left(f \delta^2\right) /\left(f+t^2\right)}\left(\frac{f}{\left(f+t^2\right)}\right)^{(f+1) / 2} \int_0^{\infty} \frac{v^f}{f !} \mathrm{e}^{-(1 / 2)\left(v-\delta t / \sqrt{f+t^2}\right)^2} \mathrm{~d} v \mathrm{~d} t\]
donde \(f\) es el número de grados de libertad (\(n-1\)), \(\delta\) es el parámetro de no centralidad (\(k \sqrt{n}\)) y \(t\) es la variable aleatoria. Esta probabilidad se basa en la función de distribución t de Student no central, se pueden encontrar más detalles en el capítulo 10 de Schilling y Neubauer (2009).
El cálculo de esta probabilidad de aceptación queda así:
Pat(p, n, k) = ccdf(NoncentralT(n-1, cquantile(Normal(),p)*sqrt(n)), k*sqrt(n))Pat (generic function with 1 method)La representación de la curva característica de operación queda así:
Pa_t = [Pat(p, n, k) for p in 0:0.005:0.2]41-element Vector{Float64}:
1.0
0.9999599292245724
0.9989007649319692
0.9938636262857532
0.9815087483518103
0.9596717321709702
0.9277552833689862
0.8864701525825365
0.8373792110645313
0.7824662021030258
0.7238066327849548
0.6633497007681655
0.6027935414546414
⋮
0.0531132263728592
0.04467527834661822
0.037491026725876786
0.03139203386176437
0.02622867333858414
0.021868911647756284
0.018196928824707026
0.01511166775538031
0.012525377493174727
0.01036219663833926
0.008556807795388699
0.007053182590451401plot(0:0.005:0.2, Pa_t, legend=false, xlabel="p", ylabel="Pₐ",
title="σ desconocida, n = $n, k = $k")Para finalizar podemos comparar las curvas caracerísticas calculadas más arriba:
plot(0:0.005:0.2, [Pa_n, Pa_t], xlabel="p", ylabel="Pₐ",
label=["σ conocida" "σ desconocida"])Dibujar las curvas características de operación de Julia no es complicado, pero requiere un cierto conocimiento de la base estadística del muestreo de aceptación y de Julia. Si no se utilizan con una cierta frecuencia es fácil olvidar como se dibujan. De hecho mi principal interés al escribir este texto es tener en un solo lugar todas estas “recetas” para poder consultarlas más adelante.