#計数規準型一回抜取検査(JIS Z 9002)の拡張
###ざっくり説明
拡張した部分は、
- JISで使われている二項分布以外に、ポアソン分布、超幾何分布を使えるようにしたこと
- JISよりも、$\alpha$と$\beta$をもうちょっと厳密に守るようにしたこと
です。
###用語など
一般的に品質保証/品質管理分野で通用する用語がわからないので、こんな意味で使ってます、という定義を書いておきます。
- $\alpha$ 生産者危険、本来合格になるはずの不適合品率のロットが不合格とされてしまう確率。
- $\beta$ 消費者危険、本来不合格となるはずの不適合品率のロットが合格とされてしまう確率。
- $AQL$ なるべく合格としたい不適合率の上限。$p_0$のこと
- $RQL$ なるべく不合格としたい不適合率の下限。 $p_1$のこと、なお$AQL<RQL$です。
- $n$ サンプルサイズ
- $c$ 合格判定個数
- $x$ 不適合品の数
###アルゴリズム
与えられた$\alpha$、$\beta$、$AQL$、$RQL$のもとで、次の式を同時に満たす$n$,$c$を探索的に求めます。
$L(AQL)$を不適合品率$AQL$のロットが合格する確率、$L(RQL)$を不適合品率$RQL$のロットが合格する確率、$P(x)$をサイズ$n$のサンプル中に不適合品が$x$個存在する確率として
\alpha = 1-L(AQL)=\sum_{x=c+1}^{n}P(x)=1-\sum_{x=0}^{c}P(x)\\
\beta = L(RQL)=\sum_{x=0}^{c}P(x)
$n$について小さい方から大きい方へ探し、最初に見つかった$n$,$c$を抜き取り方式とします。
これをjavaで書いてみました。
apache commons math3を使っています。
SamplingInspectionByAttribute.java
package QA;
import java.util.ArrayList;
import java.util.List;
import java.util.TreeMap;
import org.apache.commons.math3.distribution.AbstractIntegerDistribution;
import org.apache.commons.math3.distribution.BinomialDistribution;
import org.apache.commons.math3.distribution.HypergeometricDistribution;
import org.apache.commons.math3.distribution.PoissonDistribution;
/**
* JIS Z 9002(計数規準型抜取検査)の拡張
*
* @author gsgs
*/
public class SamplingInspectionByAttribute {
double AQL; //合格としたい品質水準
double RQL; //不合格としたい品質水準
double alpha; //生産者リスク
double beta; //消費者リスク
int max_sample; //最大サンプル数
int interval; //n探索の細かさ
DistributionType distribution;
int n; //最小サンプル数
int c; //許容基準
boolean fin = false;
/**
* 分布<br>
* BINOMIAL 二項分布 JISに同じ。通常はこれを使う<br>
* POISSON ポアソン分布 拡張。非常に大きいロットサイズかつ非常に低い不良率の時に使う<br>
* HYPERGEOMETRIC 超幾何分布 拡張。有限母集団からの非復元抽出の時に使う<br>
*/
public static enum DistributionType {
/**
* 二項分布
*/
BINOMIAL,
/**
* ポアソン分布
*/
POISSON,
/**
* 超幾何分布
*/
HYPERGEOMETRIC
}
/**
* JIS Z 9002タイプの計数規準型抜取検査<br>
* α=0.05,β=0.1,二項分布を用いて必要サンプル数(最小値)と許容範囲を求めます。
*
* @param AQL なるべく合格としたい最大不良率(0.0-1.0で与える) JISではp0
* @param RQL なるべく不合格としたい最小不良率(0.0-1.0で与える) JISではp1
*/
public SamplingInspectionByAttribute(double AQL, double RQL) {
this(AQL, RQL, 10000, DistributionType.BINOMIAL);
}
/**
* JIS Z 9002タイプの計数規準型抜取検査<br>
* α=0.05,β=0.1,指定した分布とロットサイズを用いて必要サンプル数(最小値)と許容範囲を求めます。
*
* @param AQL なるべく合格としたい最大不良率(0.0-1.0で与える) JISではp0
* @param RQL なるべく不合格としたい最小不良率(0.0-1.0で与える) JISではp1
* @param max_sample サンプルサイズ。超幾何分布の時は厳密なロットサイズを与える(他の分布の時には必要サンプル数の探索上限を与える)
* @param distribution 二項分布、ポアソン分布、超幾何分布から選ぶ
*/
public SamplingInspectionByAttribute(double AQL, double RQL, int max_sample, DistributionType distribution) {
this(AQL, RQL, 0.05, 0.1, max_sample, 1, distribution);
}
/**
* JIS Z 9002タイプの計数規準型抜取検査<br>
* 指定したパラメータ、指定した探索粒度で必要サンプル数(最小値ではないかもしれない)と許容範囲を求めます。
*
* @param AQL なるべく合格としたい最大不良率(0.0-1.0で与える) JISではp0
* @param RQL なるべく不合格としたい最小不良率(0.0-1.0で与える) JISではp1
* @param alpha 生産者危険 0.0-1.0で与える
* @param beta 消費者危険 0.0-1.0で与える
* @param max_sample サンプルサイズ。超幾何分布の時は厳密なロットサイズを与える(他の分布の時には必要サンプル数の探索上限を与える)
* @param interval
* 探索粒度。1以上かつmax_sample以下で与える。大きい数字ほど、ざっくり調べる(2以上だと必要サンプル数が最小でなくなる可能性がある、大きすぎると全数検査が必要とされやすくなる)
* @param distribution 二項分布、ポアソン分布、超幾何分布から選ぶ
*/
public SamplingInspectionByAttribute(double AQL, double RQL, double alpha, double beta, int max_sample, int interval, DistributionType distribution) {
this.AQL = AQL;
this.RQL = RQL;
this.alpha = alpha;
this.beta = beta;
this.max_sample = max_sample;
this.interval = interval;
this.interval = this.interval < 1 ? 1 : this.interval;
this.distribution = distribution;
fin = false;
for (n = interval; n < max_sample && !fin; n += interval) {
IntegerDistribution a_dist = new IntegerDistribution(max_sample, n, AQL, distribution);
IntegerDistribution b_dist = new IntegerDistribution(max_sample, n, RQL, distribution);
int min_c = a_dist.inverseCumulativeProbability(1 - alpha);
int max_c = b_dist.inverseCumulativeProbability(beta) - 1;
if (min_c <= max_c && a_dist.cumulativeProbability(min_c) >= 1 - alpha && b_dist.cumulativeProbability(min_c) <= beta) {
c = min_c;
fin = true;
break;
}
}
}
/**
* 全数検査が必要かどうかを返します
*
* @return 全数検査が必要なときtrueを返す
*/
public boolean needTotalInspection() {
return !fin;
}
/**
* 計算された最小抜取サンプル数を返す
*
* @return
*/
public int getN() {
return fin ? n : max_sample;
}
/**
* 計算された許容規準を返す。Ac
*
* @return
*/
public int getC() {
return fin ? c : (int) (AQL * max_sample);
}
/**
* サンプリングする対象番号の一例を作成する。
*
* @return サンプリングする製品の番号のリスト
*/
public Set<Integer> createSampleList() {
Set<Integer> sampleList = new TreeSet<>();
while (sampleList.size() < n) {
int x = (int) (Math.random() * max_sample);
sampleList.add(x + 1);
}
return sampleList;
}
/**
* 抜取検査計画をテキストで返す
*
* @return
*/
public String getResultReport() {
IntegerDistribution a_dist = new IntegerDistribution(max_sample, n, AQL, distribution);
IntegerDistribution b_dist = new IntegerDistribution(max_sample, n, RQL, distribution);
StringBuilder buf = new StringBuilder();
buf.append("計数抜取検査の計画\n")
.append("\n")
.append("\n")
.append("測定タイプ: 合格/不合格\n")
.append("不良率(%)によるロット品質\n")
.append("ロットサイズ: ").append(max_sample).append("\n")
.append("指定された分布: ").append(distribution.toString()).append(" Distribution\n")
.append("\n")
.append("方法\n")
.append("\n")
.append("合格品質水準 (AQL): ").append(AQL * 100).append("%\n")
.append("生産者リスク (α): ").append(alpha * 100).append("%\n")
.append("\n")
.append("不合格となる品質水準 (RQLまたはLTPD): ").append(RQL * 100).append("%\n")
.append("消費者リスク (β): ").append(beta * 100).append("%\n")
.append("\n")
.append("生成された計画\n")
.append("\n");
if (fin) {
buf.append("サンプルサイズ: ").append(n).append("\n")
.append("許容数: ").append(c).append("\n")
.append("\n")
.append("不良率(%)\t合格確率\n")
.append(AQL * 100).append("\t").append(String.format("%.4f", a_dist.cumulativeProbability(c))).append("\n")
.append(RQL * 100).append("\t").append(String.format("%.4f", b_dist.cumulativeProbability(c))).append("\n")
.append("\n");
} else {
buf.append("全数検査\n")
.append("\n");
}
return buf.toString();
}
/**
* OC曲線を返す
*
* @param min OC曲線のx軸(ロット内の不良率)の最小値を与える 0.0-1.0
* @param max OC曲線のx軸(ロット内の不良率)の最大値を与える 0.0-1.0
* @param interval OC曲線のロット内の不良率の計算間隔を与える 0.0-1.0
* @return double[OC曲線の点インデクス][x値=0,y値=1]
*/
public double[][] getOCCurve(double min, double max, double interval) {
int n = (int) ((max - min) / interval) + 1;
double[][] ret = new double[n][2];
if (needTotalInspection()) {
for (int i = 0; i < n; i++) {
ret[i][0] = min + interval * i;
ret[i][1] = ret[i][0] <= AQL ? 1 : 0;
}
} else {
for (int i = 0; i < n; i++) {
ret[i][0] = min + interval * i;
if (min + interval * i == 0d) {
ret[i][1] = 1.0d;
} else {
IntegerDistribution dist = new IntegerDistribution(max_sample, getN(), ret[i][0], distribution);
ret[i][1] = dist.cumulativeProbability(getC());
}
}
}
return ret;
}
/**
* このクラス内で使用する統計分布の抽象的なオブジェクト
*/
public static class IntegerDistribution {
AbstractIntegerDistribution dist;
public IntegerDistribution(int samples, int trial, double prob, DistributionType distribution) {
switch (distribution) {
case BINOMIAL:
dist = new BinomialDistribution(trial, prob);
break;
case POISSON:
dist = new PoissonDistribution(prob * trial);
break;
case HYPERGEOMETRIC:
dist = new HypergeometricDistribution(samples, (int) (samples * prob), trial);
break;
default:
throw new AssertionError();
}
}
public double cumulativeProbability(int x) {
return dist.cumulativeProbability(x);
}
public int inverseCumulativeProbability(double p) {
return dist.inverseCumulativeProbability(p);
}
}
}
###使い方
内部監査などで、例えば書類としてある品質記録10000枚に問題ありそうかなさそうかをまずは判断しなければいけないとき、あるベータのもとでの抜き取り検査方式をその場でリアルタイムに設計・計画できれば素早く詳細点検の必要性を判断できるのでは？と考えて、パラメタを比較的自由に設定できるように書きました。
(そのため、抜き取りNo.があれば便利だろうと、余計なおせっかいをやっています。)
動機とは異なりますが、普通の使い方での一例を次の「使い方」に書いておきます。
注意点として・・
createSampleList()を使って抜き取る試料番号を決定する場合は、max_sampleは確率分布の種類を問わず正確に与えてください。
(二項分布とポアソン分布で抜き取り検査方式が存在するかどうかだけのチェックならば、max_sampleは探索するnの上限で構いません)
使い方.java
public static void main(String[] args) {
//条件設定
double AQL = 0.005;
double RQL = 0.05;
double alpha = 0.05;
double beta = 0.1;
int max_sample = 1000;
int interval = 1;
DistributionType distribution = SamplingInspectionByAttribute.DistributionType.BINOMIAL;
//インスタンス生成
SamplingInspectionByAttribute m = new SamplingInspectionByAttribute(
AQL,
RQL,
alpha,
beta,
max_sample,
interval,
distribution
);
//結果の表示
System.out.println(m.getResultReport());
//OC曲線(本当はチャートライブラリへ突っ込んで
System.out.println("\nOC曲線");
double[][] occ = m.getOCCurve(0, RQL * 1.1, RQL / 50);
for (double[] occ1 : occ) {
System.out.printf("%.4f\t%.4f%n", occ1[0] * 100, occ1[1] * 100);
}
//便利機能：の抜き取るサンプル番号(乱数で生成、ソート済み)
System.out.println("\n抜き取りNo\n");
Set<Integer> sampleList = m.createSampleList();
int cnt = 0;
for (Integer samp : sampleList) {
System.out.print("□ " + samp + "\t");
cnt++;
if (cnt % 5 == 0) {
System.out.println("");
}
}
System.out.println("------");
}
###実行結果
run:
計数抜取検査の計画
測定タイプ: 合格/不合格
不良率(%)によるロット品質
ロットサイズ: 1000
指定された分布: BINOMIAL Distribution
方法
合格品質水準 (AQL): 0.5%
生産者リスク (α): 5.0%
不合格となる品質水準 (RQLまたはLTPD): 5.0%
消費者リスク (β): 10.0%
生成された計画
サンプルサイズ: 105
許容数: 2
不良率(%) 合格確率
0.5 0.9839
5.0 0.0992
OC曲線
0.0000 100.0000
0.1000 99.9826
0.2000 99.8712
0.3000 99.5970
0.4000 99.1141
0.5000 98.3947
0.6000 97.4253
0.7000 96.2039
0.8000 94.7366
0.9000 93.0363
1.0000 91.1201
1.1000 89.0084
1.2000 86.7235
1.3000 84.2887
1.4000 81.7276
1.5000 79.0632
1.6000 76.3182
1.7000 73.5137
1.8000 70.6698
1.9000 67.8051
2.0000 64.9366
2.1000 62.0796
2.2000 59.2480
2.3000 56.4540
2.4000 53.7081
2.5000 51.0198
2.6000 48.3967
2.7000 45.8454
2.8000 43.3714
2.9000 40.9789
3.0000 38.6710
3.1000 36.4503
3.2000 34.3182
3.3000 32.2756
3.4000 30.3225
3.5000 28.4587
3.6000 26.6832
3.7000 24.9946
3.8000 23.3914
3.9000 21.8715
4.0000 20.4328
4.1000 19.0727
4.2000 17.7887
4.3000 16.5781
4.4000 15.4381
4.5000 14.3657
4.6000 13.3582
4.7000 12.4126
4.8000 11.5261
4.9000 10.6957
5.0000 9.9187
5.1000 9.1923
5.2000 8.5139
5.3000 7.8807
5.4000 7.2904
5.5000 6.7404
抜き取りNo
□ 2 □ 8 □ 12 □ 40 □ 45
□ 57 □ 58 □ 66 □ 74 □ 85
□ 87 □ 96 □ 107 □ 113 □ 120
□ 122 □ 134 □ 142 □ 146 □ 150
□ 163 □ 172 □ 175 □ 186 □ 205
□ 207 □ 213 □ 237 □ 247 □ 251
□ 253 □ 261 □ 276 □ 278 □ 286
□ 293 □ 306 □ 307 □ 310 □ 316
□ 319 □ 320 □ 337 □ 366 □ 374
□ 385 □ 392 □ 407 □ 419 □ 421
□ 429 □ 457 □ 471 □ 479 □ 518
□ 528 □ 533 □ 546 □ 564 □ 572
□ 586 □ 588 □ 598 □ 599 □ 611
□ 617 □ 628 □ 629 □ 649 □ 650
□ 674 □ 677 □ 683 □ 689 □ 691
□ 692 □ 710 □ 713 □ 722 □ 739
□ 744 □ 752 □ 756 □ 761 □ 794
□ 807 □ 809 □ 830 □ 843 □ 863
□ 881 □ 888 □ 889 □ 898 □ 902
□ 906 □ 913 □ 917 □ 931 □ 949
□ 964 □ 984 □ 993 □ 995 □ 1000
------
ビルド成功(合計時間: 0秒)
###お願い
実用とするにはGUIが必要ですので、お好きにかぶせてください。
間違いに気づかれた方は、どうか、教えてください。
n,cを探索せずとも、一発で求められるよ、という方もぜひ教えてください。