diff --git a/src/main/java/info/istlab/Zemi01/miuramath/Utest.java b/src/main/java/info/istlab/Zemi01/miuramath/Utest.java
new file mode 100644
index 0000000..724ee3b
--- /dev/null
+++ b/src/main/java/info/istlab/Zemi01/miuramath/Utest.java
@@ -0,0 +1,437 @@
+package info.istlab.Zemi01.miuramath;
+
+import java.util.Arrays;
+
+import org.apache.commons.math3.distribution.TDistribution;
+
+public class Utest {
+    public static void main(String[] args) {
+        double[] d1 = { 4, 3, 3, 3, 2, 3, 4, 5, 4, 3 };
+        double[] d2 = { 5, 4, 5, 4, 4, 5, 4, 3, 5, 4 };
+
+        System.out.println(Arrays.toString(utest(d1, d2)));
+    }
+
+    public static double[] utest(double[] d1, double[] d2) {
+        boolean isDebug = false;
+        double[] sorted = merged_and_sorted(d1, d2);
+        if (isDebug) System.out.println(Arrays.toString(sorted));
+
+        double[] rnkd1 = ranks(d1, sorted); // rank配列
+        double[] rnkd2 = ranks(d2, sorted);
+        if (isDebug) System.out.println(Arrays.toString(rnkd1));
+        if (isDebug) System.out.println(Arrays.toString(rnkd2));
+
+        double[] samed1 = samenum(d1, sorted); // 同位の数の配列
+        double[] samed2 = samenum(d2, sorted); // 同位の数の配列
+        if (isDebug) System.out.println(Arrays.toString(samed1));
+        if (isDebug) System.out.println(Arrays.toString(samed2));
+
+        double[] fo1 = fixedorder(rnkd1, samed1); // 調整済み
+        double[] fo2 = fixedorder(rnkd2, samed2); // 調整済み
+        if (isDebug) System.out.println(Arrays.toString(fo1));
+        if (isDebug) System.out.println(Arrays.toString(fo2));
+
+        // (同位数^3-同位数)/同位数
+        double[] ddd1 = sameordpow(samed1);
+        double[] ddd2 = sameordpow(samed2);
+        if (isDebug) System.out.println(Arrays.toString(ddd1));
+        if (isDebug) System.out.println(Arrays.toString(ddd2));
+
+        double K18 = sum(ddd1) + sum(ddd2);
+        int B18 = d1.length;
+        int C18 = d2.length;
+        int B19 = B18 + C18;
+        double H18 = sum(fo1);
+        double I18 = sum(fo2);
+
+        double tempU1 = B18 * C18 + (B18 * (B18 + 1) / 2) - H18;
+        double tempU2 = B18 * C18 + (B18 * (B18 + 1) / 2) - I18;
+        double B21 = (tempU1 < tempU2) ? tempU1 : tempU2; // U
+        if (isDebug) System.out.println("U = "+B21);
+        double B22 = (tempU1 + tempU2) / 2;// aveU
+
+        double B23 = B18 * C18 * (Math.pow(B19, 3) - B19 - K18) / (12 * (Math.pow(B19, 2) - B19));
+        double B24 = Math.abs(B21 - B22) / Math.sqrt(B23);
+
+        var nomdist = new org.apache.commons.math3.distribution.NormalDistribution();
+
+        double B25 = (1 - nomdist.cumulativeProbability(B24)) * 2;
+        return new double[]{B21, B25}; //U, p
+    }
+
+    public static double[] merged_and_sorted(double[] d1, double[] d2) {
+        int n1 = d1.length;
+        int n2 = d2.length;
+        double[] merged = new double[n1 + n2];
+        for (int i = 0; i < n1; i++) {
+            merged[i] = d1[i];
+        }
+        for (int i = 0; i < n2; i++) {
+            merged[i + n1] = d2[i];
+        }
+        Arrays.sort(merged); // ここでソート。ただし昇順
+        return merged;
+    }
+
+    /**
+     * refと同じサイズの、ランク配列を返す。
+     * 
+     * @param ref       元データ
+     * @param incsorted 昇順ソート済みの配列
+     * @return
+     */
+    public static double[] ranks(double[] ref, double[] incsorted) {
+        double[] rnk = new double[ref.length];
+        int isn = incsorted.length;
+
+        for (int i = 0; i < ref.length; i++) {
+            // incsortedの後ろから、ref[i]が何個目かを探す
+            for (int j = 1; j <= incsorted.length; j++) {
+                if (incsorted[isn - j] == ref[i]) {
+                    rnk[i] = j;
+                    break;
+                }
+            }
+        }
+        return rnk;
+    }
+
+    /**
+     * refと同じサイズの、同位の数の配列を返す。
+     * 
+     * @param ref       元データ
+     * @param incsorted 昇順ソート済みの配列
+     * @return
+     */
+    public static double[] samenum(double[] ref, double[] incsorted) {
+        double[] same = new double[ref.length];
+        int isn = incsorted.length;
+
+        for (int i = 0; i < ref.length; i++) {
+            // incsortedの後ろから探して、ref[i]と同じなら、カウントアップ
+            for (int j = 1; j <= incsorted.length; j++) {
+                if (incsorted[isn - j] == ref[i]) {
+                    same[i]++;
+                }
+            }
+        }
+        return same;
+    }
+
+    /**
+     * 同位調整済み順位の配列を返す
+     * 
+     * @param rankd
+     * @param samed
+     * @return
+     */
+    public static double[] fixedorder(double[] rankd, double[] samed) {
+        double[] fixedorder = new double[rankd.length];
+        for (int i = 0; i < rankd.length; i++) {
+            double rnk = rankd[i];
+            double sam = samed[i];
+            fixedorder[i] = (sam * (rnk + rnk + sam - 1) / 2) / sam;
+        }
+        return fixedorder;
+    }
+
+    /**
+     * (同位数^3-同位数)/同位数
+     * 
+     * @param d
+     * @return
+     */
+    public static double[] sameordpow(double[] d) {
+        double[] ret = new double[d.length];
+        for (int i = 0; i < d.length; i++) {
+            double x = d[i];
+            ret[i] = (Math.pow(x, 3) - x) / x;
+        }
+        return ret;
+    }
+
+    /**
+     * 合計
+     */
+    public static double sum(double[] ary) {
+        double tmp = 0;
+        for (int i = 0; i < ary.length; i++) {
+            tmp += ary[i];
+        }
+        return tmp;
+    }
+
+    /**
+     * 平均
+     */
+    public static double avg(double[] ary) {
+        double sum = sum(ary);
+        return sum / ary.length;
+    }
+
+    /**
+     * 偏差：平均値との差の2乗
+     */
+    public static double hensa(double[] ary) {
+        double avg = avg(ary);
+        double tmp = 0;
+        for (int i = 0; i < ary.length; i++) {
+            tmp += (ary[i] - avg) * (ary[i] - avg);
+        }
+        return tmp;
+    }
+
+    /**
+     * 不偏分散： 偏差を、データ数n-1 で割ったもの vars
+     */
+    public static double var(double[] ary) {
+        double hensa = hensa(ary);
+        return hensa / (ary.length - 1);
+    }
+
+    /**
+     * 標本分散： 偏差を、データ数n で割ったもの
+     */
+    public static double varp(double[] ary) {
+        double hensa = hensa(ary);
+        return hensa / (ary.length);
+    }
+
+    /**
+     * 差の配列
+     * 
+     * @param ary1
+     * @param ary2
+     * @return
+     */
+    public static double[] diffary(double[] ary1, double[] ary2) {
+        double[] ret = new double[ary1.length];
+        for (int i = 0; i < ary1.length; i++) {
+            ret[i] = ary1[i] - ary2[i];
+        }
+        return ret;
+    }
+
+    /**
+     * 対応のあるt検定の t値
+     */
+    public static double paired_t(double[] ary1, double[] ary2) {
+        // 差の平均値 / sqrt(差の不偏分散 / データ数)
+        double avg_sa = avg(diffary(ary1, ary2));
+        double var_sa = var(diffary(ary1, ary2));
+        return avg_sa / Math.sqrt(var_sa / ary1.length);
+    }
+
+    /**
+     * 対応がない(なくてもよい) t検定 Welchのt
+     */
+    public static double welch_t(double[] ary1, double[] ary2) {
+        // 不偏分散/データ数
+        double stddev1 = var(ary1) / ary1.length;
+        double stddev2 = var(ary2) / ary2.length;
+        double diff_avg = avg(ary1) - avg(ary2);
+        return diff_avg / Math.sqrt(stddev1 + stddev2);
+    }
+
+    public static double welch_degree_of_freedom(double[] ary1, double[] ary2) {
+        // 不偏分散/データ数
+        double stddev1 = var(ary1) / ary1.length;
+        double stddev2 = var(ary2) / ary2.length;
+        double a1 = Math.pow(stddev1 + stddev2, 2);
+        double b1 = Math.pow(stddev1, 2) / (ary1.length - 1) + Math.pow(stddev2, 2) / (ary2.length - 1);
+        return a1 / b1;
+    }
+
+    /**
+     * TDIST by commons-math3 TDistribution
+     * 
+     * @param x
+     * @param df
+     * @param tails
+     * @return
+     */
+    public static double tdist(double x, double deg_free, int tails) {
+        TDistribution tdistobj = new TDistribution(deg_free);
+        double cdf = tdistobj.cumulativeProbability(x);
+        double ret = (1 - cdf) * tails;
+        // //test
+        // double test = tdist_old(x, deg_free, tails);
+        // if (test != ret) System.out.println("DIFF ret = "+ret+" test "+test );
+        return ret;
+    }
+
+    /**
+     * 
+     * @param x
+     * @param deg_free
+     * @param tails
+     * @return
+     */
+    public static double tdist_old(double x, double deg_free, int tails) {
+        return (1 - studentsCdf(x, deg_free)) * tails;
+    }
+
+    /**
+     * 
+     * http://www.java2s.com/example/java-api/java/lang/math/exp-1-32.html
+     * 
+     * 
+     * Calculates the probability from -INF to X under Student's Distribution
+     * Ported to PHP from Javascript implementation found at
+     * http://www.math.ucla.edu/~tom/distributions/tDist.html
+     * // w w w .j av a 2 s. c o m
+     * 
+     * @param x
+     * @param deg_free
+     * @return
+     */
+    public static double studentsCdf(double x, double deg_free) {
+        if (deg_free <= 0) {
+            throw new IllegalArgumentException("The degrees of freedom need to be positive.");
+        }
+
+        double A = deg_free / 2.0;
+        double S = A + 0.5;
+        double Z = deg_free / (deg_free + x * x);
+        double BT = Math.exp(logGamma(S) - logGamma(0.5) - logGamma(A) + A * Math.log(Z) + 0.5 * Math.log(1.0 - Z));
+        double betacdf;
+        if (Z < (A + 1.0) / (S + 2.0)) {
+            betacdf = BT * betinc(Z, A, 0.5);
+        } else {
+            betacdf = 1 - BT * betinc(1.0 - Z, 0.5, A);
+        }
+
+        double tcdf;
+        if (x < 0) {
+            tcdf = betacdf / 2.0;
+        } else {
+            tcdf = 1.0 - betacdf / 2.0;
+        }
+        return tcdf;
+    }
+
+    /**
+     * Log Gamma Function
+     * 
+     * @param Z
+     * @return
+     */
+    public static double logGamma(double Z) {
+        double S = 1.0 + 76.18009173 / Z - 86.50532033 / (Z + 1.0) + 24.01409822 / (Z + 2.0)
+                - 1.231739516 / (Z + 3.0) + 0.00120858003 / (Z + 4.0) - 0.00000536382 / (Z + 5.0);
+        double LG = (Z - 0.5) * Math.log(Z + 4.5) - (Z + 4.5) + Math.log(S * 2.50662827465);
+
+        return LG;
+    }
+
+    /**
+     * Internal function used by StudentCdf
+     * 
+     * @param x
+     * @param A
+     * @param B
+     * @return
+     */
+    protected static double betinc(double x, double A, double B) {
+        double A0 = 0.0;
+        double B0 = 1.0;
+        double A1 = 1.0;
+        double B1 = 1.0;
+        double M9 = 0.0;
+        double A2 = 0.0;
+        while (Math.abs((A1 - A2) / A1) > 0.00000001) {
+            A2 = A1;
+            double C9 = -(A + M9) * (A + B + M9) * x / (A + 2.0 * M9) / (A + 2.0 * M9 + 1.0);
+            A0 = A1 + C9 * A0;
+            B0 = B1 + C9 * B0;
+            M9 = M9 + 1;
+            C9 = M9 * (B - M9) * x / (A + 2.0 * M9 - 1.0) / (A + 2.0 * M9);
+            A1 = A0 + C9 * A1;
+            B1 = B0 + C9 * B1;
+            A0 = A0 / B1;
+            B0 = B0 / B1;
+            A1 = A1 / B1;
+            B1 = 1.0;
+        }
+        return A1 / A;
+    }
+
+    /**
+     * うまくいかないので使わない
+     * 
+     * @param x
+     * @param v
+     * @return
+     */
+    public static double studentsCdf2(double x, double v) {
+        double xL = x;
+        double xR = 100000;
+        int n = 10000;
+        double goukei = 0.0;
+        double a, b, xc, tmp;
+        double dx = (xR - xL) / n;
+        double B = simpson_betafunc(0.5, v / 2);
+        // double B = Beta.regularizedBeta(0.5, v / 2);
+        for (int i = 0; i < n; i++) {
+            a = xL + i * dx;
+            b = a + dx;
+            xc = a + dx / 2.0;
+            tmp = 0.0;
+            tmp = tmp + densityfunc(a, v, B) / 3;
+            tmp = tmp + densityfunc(xc, v, B) * 4 / 3;
+            tmp = tmp + densityfunc(b, v, B) / 3;
+            tmp = tmp * dx / 2;
+            goukei += tmp;
+        }
+        return goukei;
+    }
+
+    private static double densityfunc(double x, double v, double B) {
+        return Math.pow(1.0 + (x * x / v), -(v + 1) / 2) / (Math.sqrt(v) * B);
+    }
+
+    /**
+     * うまくいかないので使わない
+     * 
+     * @param u
+     * @param v
+     * @return
+     */
+    public static double simpson_betafunc(double u, double v) {
+        double xL = 0.0;
+        double xR = 1.0;
+        int n = 100;
+        double goukei = 0.0;
+        double a, b, xc, tmp;
+        double dx = (xR - xL) / n;
+        for (int i = 0; i < n; i++) {
+            a = xL + i * dx;
+            b = a + dx;
+            xc = a + dx / 2.0;
+            tmp = 0.0;
+            tmp = tmp + bfunc(a, u, v) / 3;
+            tmp = tmp + bfunc(xc, u, v) * 4 / 3;
+            tmp = tmp + bfunc(b, u, v) / 3;
+            tmp = tmp * dx / 2;
+            goukei += tmp;
+            System.out.println(i + " " + tmp + " " + goukei);
+        }
+        return goukei;
+    }
+
+    /**
+     * うまくいかないので使わない
+     * 
+     * @param x
+     * @param u
+     * @param v
+     * @return
+     */
+    private static double bfunc(double x, double u, double v) {
+        double d = Math.pow(x, u - 1) * Math.pow(1.0 - x, v - 1);
+        System.out.println("bfunc " + x + " " + u + " " + v + " = " + d);
+        return d;
+    }
+
+}