GC補正
やっつけ仕事でGC補正をしたいとき。
public class gC { public double lambdaori; public double lambda; public double[] correctedPs; public static double chimedian = 0.454936425; public double minpOri; public double minpCorr; public double minpCorrplus; public gC(double[] ps){ correctedPs = StatUtils.MiscUtil.DeepCopyDouble1(ps); double[] tmp = StatUtils.MiscUtil.DeepCopyDouble1(ps); //System.out.println("presort"); StatUtils.MiscUtil.quickSortDouble(tmp,0,tmp.length-1); //System.out.println("postsort"); minpOri=tmp[0]; //System.out.println("minp="+minpOri); double median = StatUtils.MiscUtil.quantileFromSortedList2(tmp, 0.5); //System.out.println("median="+median); org.apache.commons.math.distribution.ChiSquaredDistributionImpl ch = new org.apache.commons.math.distribution.ChiSquaredDistributionImpl(1); try{ double medchi = ch.inverseCumulativeProbability(1-median); lambdaori = medchi/chimedian; if(lambdaori<1){ lambda=1; }else{ lambda=lambdaori; } //System.out.println("lambdaori=\t"+lambdaori+"\tlambda\t"+lambda); minpCorr=1-ch.cumulativeProbability(ch.inverseCumulativeProbability(1-minpOri)/lambda); minpCorrplus=1-ch.cumulativeProbability(ch.inverseCumulativeProbability(1-minpOri)/lambdaori); //System.out.println("minpOr\t"+minpOri+"\tminpCorr\t"+minpCorr); } catch (org.apache.commons.math.MathException ex) { } } public gC(double[] ps,String corrop){ correctedPs = StatUtils.MiscUtil.DeepCopyDouble1(ps); double[] tmp = StatUtils.MiscUtil.DeepCopyDouble1(ps); //System.out.println("presort"); StatUtils.MiscUtil.quickSortDouble(tmp,0,tmp.length-1); //System.out.println("postsort"); minpOri=tmp[0]; //System.out.println("minp="+minpOri); double median = StatUtils.MiscUtil.quantileFromSortedList2(tmp, 0.5); //System.out.println("median="+median); org.apache.commons.math.distribution.ChiSquaredDistributionImpl ch = new org.apache.commons.math.distribution.ChiSquaredDistributionImpl(1); try{ double medchi = ch.inverseCumulativeProbability(1-median); lambdaori = medchi/chimedian; if(lambdaori<1){ lambda=1; }else{ lambda=lambdaori; } //System.out.println("lambdaori=\t"+lambdaori+"\tlambda\t"+lambda); minpCorr=1-ch.cumulativeProbability(ch.inverseCumulativeProbability(1-minpOri)/lambda); minpCorrplus=1-ch.cumulativeProbability(ch.inverseCumulativeProbability(1-minpOri)/lambdaori); if(corrop.equals("c")){ for(int i=0;i<correctedPs.length;i++){ correctedPs[i]=1-ch.cumulativeProbability(ch.inverseCumulativeProbability(1-ps[i])/lambda); } } //System.out.println("minpOr\t"+minpOri+"\tminpCorr\t"+minpCorr); } catch (org.apache.commons.math.MathException ex) { } } public gC(double[] tmp,int type,org.apache.commons.math.distribution.ChiSquaredDistributionImpl ch){ if(type==9){ /* * 入力Pを保存しない */ } //correctedPs = StatUtils.MiscUtil.DeepCopyDouble1(ps); //double[] tmp = StatUtils.MiscUtil.DeepCopyDouble1(ps); //System.out.println("presort"); StatUtils.MiscUtil.quickSortDouble(tmp,0,tmp.length-1); //System.out.println("postsort"); minpOri=tmp[0]; //System.out.println("minp="+minpOri); double median = StatUtils.MiscUtil.quantileFromSortedList2(tmp, 0.5); //System.out.println("median="+median); //org.apache.commons.math.distribution.ChiSquaredDistributionImpl ch = //new org.apache.commons.math.distribution.ChiSquaredDistributionImpl(1); try{ double medchi = ch.inverseCumulativeProbability(1-median); lambdaori = medchi/chimedian; if(lambdaori<1){ lambda=1; }else{ lambda=lambdaori; } //System.out.println("lambdaori=\t"+lambdaori+"\tlambda\t"+lambda); double tmp1 = ch.inverseCumulativeProbability(1-minpOri)/lambda; double tmp2 = ch.inverseCumulativeProbability(1-minpOri)/lambdaori; //System.out.println("tmp1="+tmp1+" tmp2="+tmp2); minpCorr=1-ch.cumulativeProbability(tmp1); minpCorrplus=1-ch.cumulativeProbability(tmp2); //System.out.println("minpOr\t"+minpOri+"\tminpCorr\t"+minpCorr); } catch (org.apache.commons.math.MathException ex) { } } public gC(double[] ps,int a, int b){ double[] partPs = new double[b-a+1]; for(int i=a;i<=b;i++){ partPs[i-a]=ps[i]; } gC tmpgC = new gC(partPs); correctedPs =tmpgC.correctedPs; minpOri=tmpgC.minpOri; minpCorr=tmpgC.minpCorr; lambda=tmpgC.lambda; lambdaori=tmpgC.lambdaori; minpCorrplus=tmpgC.minpCorrplus; } public gC(double[] ps,int a, int b,int nosave,org.apache.commons.math.distribution.ChiSquaredDistributionImpl ch){ double[] partPs = new double[b-a+1]; for(int i=a;i<=b;i++){ partPs[i-a]=ps[i]; } gC tmpgC = new gC(partPs,nosave,ch); correctedPs =tmpgC.correctedPs; minpOri=tmpgC.minpOri; minpCorr=tmpgC.minpCorr; lambda=tmpgC.lambda; lambdaori=tmpgC.lambdaori; minpCorrplus=tmpgC.minpCorrplus; } public String StringgC(gC g,String sep1,String sep2){ String ret =""; ret += "lambdaori"+sep1+g.lambdaori+sep1; ret += "lambda"+sep1+g.lambda+sep1; ret += "minpOri"+sep1+g.minpOri+sep1; ret += "minpCorr"+sep1+g.minpCorr+sep1; ret += "minpCorrplus"+sep1+g.minpCorrplus+sep1; return ret; } public void PrintgC(gC g,String sep1,String sep2){ String ret = StringgC(g,sep1,sep2); System.out.println(ret); } public void PrintCorrectedPs(String sep1,String sep2){ for(int i=0;i< correctedPs.length;i++){ System.out.print(correctedPs[i]+sep2); } } }