2006-10-14

カイ自乗統計量

独立性検定 Java

複数のカテゴリカル変数からなる多次元分割表を考える。k変数あり、それぞれがdk個のカテゴリを持つとき、分割表のセル数は $￥prod dk$ である。このセルについて、独立仮説に基づいた期待値を求め、観測値と期待値とから算出されるカイ自乗統計量は、多次元分割表の独立性検定を行うときの、適当な統計量である。たとえば、パーミュテーションテストを行う場合など。



public static void main(String[] args) {

		int[] data = {2,2,4,4,2,2,4,4};

		int[] dim = {2,2,2};

		int[] data2 = {1,2,3,4};

		int[] dim2 = {2,2};

		int[] data3 = {2,2,2,4,4,4};

		int[] dim3 = {3,2};

		double multChi = chiMultDim(data,dim);

		System.out.println("multChi " + multChi);

		double multChi2 = chiMultDim(data2,dim2);

		System.out.println("multChi2 " + multChi2);

		

		double multChi3 = chiMultDim(data3,dim3);

		System.out.println("multChi3 " + multChi3);

	}

public static double chiMultDim(int[] d, int[] dim){

		double ret =0;

		

		int[][] lastMarginal = LastMarginal(d, dim);

		int sum = Calculator.sumInt(d); 

		

		double[] exp = new double[d.length];

		for(int i=0;i<exp.length;i++){

			exp[i] =1;

			for(int j=0;j<lastMarginal[i].length;j++){

				exp[i] *= (double)lastMarginal[i][j]/(double)sum;

			}

			exp[i] *= (double)sum;

		}

		for(int i=0;i<d.length;i++){

			ret +=((double)d[i]-exp[i])*((double)d[i]-exp[i])/exp[i];

		}

		return ret;

		

		

	}

public static int[][] LastMarginal(int[] d, int[] dim){

		int[][] ret = new int[d.length][2];//min, max

		//int[][] marginal = new int[dim.length][0];

		int[][] marginal = new int[d.length][dim.length];

		//for(int i=0;i<marginal.length;i++){

			//marginal[i] = new int[dim[i]];

		//}

		int[][] max = new int[d.length][dim.length];

		int[] keta = new int[dim.length];

		int[] keta2 = new int[dim.length];

		for(int i=0;i<keta.length;i++){

			keta[i]=0;

			keta2[i]=0;

		}

		for(int i=0;i<d.length;i++){

			for(int k=0;k<keta.length;k++){

				//keta[i]=0;

				keta2[k]=0;

			}

			for(int j=0;j<marginal.length;j++){

				

				for(int k=0;k<marginal[j].length;k++){

					//System.out.println("keta2 keta k" + + k +" " + keta2[k] + " " + keta[k]);

					

					if(keta2[k]==keta[k]){

						//System.out.println("YES keta2 keta " + keta2[k] + " " + keta[k]);

						

						marginal[j][k]+=d[i];

					}

				}

				keta2[0]++;

				boolean success = false;

				for(int m=0;m<keta2.length-1;m++){

					if(keta2[m]==dim[m]){

						keta2[m]=0;

						

						keta2[m+1]++;

						

						

					}else{

						success = true;

						break;

					}

				}

				if(keta2[keta2.length-1]<dim[keta2.length-1]){

					success = true;

				}

			}

			

			/*

			for(int j=0;j<keta.length;j++){

				marginal[j][keta[j]]+=d[i];

			}

			*/

			

			keta[0]++;

			boolean success = false;

			for(int j=0;j<keta.length-1;j++){

				if(keta[j]==dim[j]){

					keta[j]=0;

					

					keta[j+1]++;

					

					

				}else{

					success = true;

					break;

				}

			}

			if(keta[keta.length-1]<dim[keta.length-1]){

				success = true;

			}

			

		}

		String Stmarginal ="";

		for(int i=0;i<marginal.length;i++){

			for(int j=0;j<marginal[i].length;j++){

				Stmarginal += marginal[i][j] +"\t";

			}

			Stmarginal += "\n";

		}

		System.out.println(Stmarginal);

		return marginal;

		/*

		for(int i=0;i<max.length;i++){

			ret[i][0]=0;

			ret[i][1]=max[i][0];

		}

		for(int i=0;i<keta.length;i++){

			keta[i]=0;

		}

		for(int i=0;i<d.length;i++){

			ret[i][1]=marginal[i][0];

			for(int j=1;j<marginal[i].length;j++){

				if(ret[i][1]>marginal[i][j]){

					ret[i][1] = marginal[i][j];

				}

			}

			

			keta[0]++;

			boolean success = false;

			for(int j=0;j<keta.length-1;j++){

				if(keta[j]==dim[j]){

					keta[j]=0;

					

					keta[j+1]++;

					

					

				}else{

					success = true;

					break;

				}

			}

			if(keta[keta.length-1]<dim[keta.length-1]){

				success = true;

			}

			

		}

		
	

		return ret;

		*/

		

	}