[BACK]Return to gtt_ekn-ja.texi CVS log [TXT][DIR] Up to [local] / OpenXM / src / asir-contrib / packages / doc / gtt_ekn

Annotation of OpenXM/src/asir-contrib/packages/doc/gtt_ekn/gtt_ekn-ja.texi, Revision 1.6

1.6     ! takayama    1: %% $OpenXM: OpenXM/src/asir-contrib/packages/doc/gtt_ekn/gtt_ekn-ja.texi,v 1.5 2017/03/03 12:02:32 takayama Exp $
        !             2: %% ptex -kanji=utf8 gtt_ekn.texi   (.texi までつける. platex でなく ptex)
        !             3: %% 以下コメントは @comment で始める.  \input texinfo 以降は普通の tex 命令は使えない.
1.1       takayama    4: \input texinfo
                      5: @iftex
                      6: @catcode`@#=6
                      7: @def@fref#1{@xrefX[#1,,@code{#1},,,]}
                      8: @def@b#1{{@bf@gt #1}}
                      9: @catcode`@#=@other
                     10: @end iftex
                     11: @overfullrule=0pt
                     12: @c -*-texinfo-*-
                     13: @comment %**start of header
1.6     ! takayama   14: @comment --- おまじない終り ---
1.1       takayama   15:
1.6     ! takayama   16: @comment --- GNU info ファイルの名前 ---
1.1       takayama   17: @setfilename xyzman
                     18:
1.6     ! takayama   19: @comment --- タイトル ---
        !            20: @settitle 2元分割表HGM
1.1       takayama   21:
                     22: @comment %**end of header
                     23: @comment %@setchapternewpage odd
                     24:
1.6     ! takayama   25: @comment --- おまじない ---
1.1       takayama   26: @ifinfo
                     27: @macro fref{name}
                     28: @ref{\name\,,@code{\name\}}
                     29: @end macro
                     30: @end ifinfo
                     31:
                     32: @iftex
                     33: @comment @finalout
                     34: @end iftex
                     35:
                     36: @titlepage
1.6     ! takayama   37: @comment --- おまじない終り ---
1.1       takayama   38:
1.6     ! takayama   39: @comment --- タイトル, バージョン, 著者名, 著作権表示 ---
        !            40: @title 2元分割表HGM関数
        !            41: @subtitle Risa/Asir 2元分割表HGM関数説明書
        !            42: @subtitle 1.1 版
        !            43: @subtitle 2017 年 3 月 3 日
1.1       takayama   44:
                     45: @author  by Y.Goto, Y.Tachibana, N.Takayama
                     46: @page
                     47: @vskip 0pt plus 1filll
                     48: Copyright @copyright{} Risa/Asir committers
                     49: 2004--2010. All rights reserved.
                     50: @end titlepage
                     51:
1.6     ! takayama   52: @comment --- おまじない ---
1.1       takayama   53: @synindex vr fn
1.6     ! takayama   54: @comment --- おまじない終り ---
1.1       takayama   55:
1.6     ! takayama   56: @comment --- @node は GNU info, HTML 用 ---
        !            57: @comment --- @node  の引数は node-name,  next,  previous,  up ---
1.1       takayama   58: @node Top,, (dir), (dir)
                     59:
1.6     ! takayama   60: @comment --- @menu は GNU info, HTML 用 ---
        !            61: @comment --- chapter 名を正確に並べる ---
        !            62: @comment --- この文書では chapter XYZ, Chapter Index がある.
        !            63: @comment ---  Chapter XYZ には section XYZについて, section XYZに関する関数がある.
1.1       takayama   64: @menu
1.6     ! takayama   65: * 2元分割表HGMの関数説明書について::
        !            66: * 2元分割表HGMの関数::
        !            67: * modular計算
1.1       takayama   68: * Index::
                     69: @end menu
                     70:
1.6     ! takayama   71: @comment --- chapter の開始 ---
        !            72: @comment --- 親 chapter 名を正確に. 親がない場合は Top ---
        !            73: @node 2元分割表HGMの関数説明書について,,, Top
        !            74: @chapter 2元分割表HGMの関数説明書について
        !            75:
        !            76: この説明書では
        !            77: HGM(holonomic gradient method) を用いた2元分割表の関数について説明する.
        !            78: ChangeLog の項目は www.openxm.org の cvsweb で
        !            79: ソースコードを読む時の助けになる情報が書かれている.
1.1       takayama   80:
1.6     ! takayama   81: 本文中で引用している文献を列挙する.
1.1       takayama   82: @itemize @bullet
                     83: @item [GM2016]
                     84: Y.Goto, K.Matsumoto, Pfaffian equations and contiguity relations of the hypergeometric function of type (k+1,k+n+2) and their applications, arxiv:1602.01637 (version 1)
                     85: @item [T2016]
1.6     ! takayama   86: Y.Tachibana, 差分ホロノミック勾配法のモジュラーメソッドによる計算の高速化,
        !            87: 2016, 神戸大学修士論文.
1.1       takayama   88: @item [GTT2016]
1.6     ! takayama   89: Y.Goto, Y.Tachibana, N.Takayama, 2元分割表に対する差分ホロノミック勾配法の実装,
        !            90: 数理研講究録(掲載予定).
1.1       takayama   91: @item [TKT2015]
                     92: N.Takayama, S.Kuriki, A.Takemura,
                     93:          $A$-hypergeometric distributions and Newton polytopes.
                     94:          arxiv:1510.02269
                     95: @end itemize
                     96:
1.6     ! takayama   97: このマニュアルで説明する関数を用いたプログラム例は
1.1       takayama   98: gtt_ekn/test-t1.rr
1.6     ! takayama   99: など.
1.1       takayama  100:
1.4       takayama  101:
1.6     ! takayama  102: @node 2元分割表HGMの関数,,, Top
        !           103: @chapter 2元分割表HGMの関数
1.1       takayama  104:
1.6     ! takayama  105: @comment --- section ``実験的関数'' の subsection xyz_abc
        !           106: @comment --- subsection xyz_pqr xyz_stu がある.
1.1       takayama  107: @menu
                    108: * gtt_ekn.gmvector::
                    109: * gtt_ekn.nc::
                    110: * gtt_ekn.lognc::
                    111: * gtt_ekn.expectation::
                    112: * gtt_ekn.setup::
                    113: * gtt_ekn.upAlpha::
1.5       takayama  114: * gtt_ekn.cmle::
1.1       takayama  115: @end menu
                    116:
1.6     ! takayama  117: @node 超幾何関数E(k,n),,, 2元分割表HGMの関数
        !           118: @section 超幾何関数E(k,n)
1.1       takayama  119:
                    120: @comment **********************************************************
1.6     ! takayama  121: @comment --- ◯◯◯◯  の説明
        !           122: @comment --- 個々の関数の説明の開始 ---
        !           123: @comment --- section 名を正確に ---
        !           124: @node gtt_ekn.gmvector,,, 超幾何関数E(k,n)
1.1       takayama  125: @subsection @code{gtt_ekn.gmvector}
1.6     ! takayama  126: @comment --- 索引用キーワード
1.1       takayama  127: @findex gtt_ekn.gmvector
                    128:
                    129: @table @t
                    130: @item gtt_ekn.gmvector(@var{beta},@var{p})
1.6     ! takayama  131: :: 周辺和 @var{beta}, セルの確率 @var{p} の二元分割表に付随する超幾何関数
        !           132: E(k,n) の値およびその微分の値を戻す.
1.1       takayama  133: @item gtt_ekn.ekn_cBasis_2(@var{beta},@var{p})
1.6     ! takayama  134: の別名である.
1.1       takayama  135: @end table
                    136:
1.6     ! takayama  137: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  138: @table @var
                    139: @item return
1.6     ! takayama  140: ベクトル, 超幾何関数の値とその微分. 詳しくは下記.
1.1       takayama  141: @item beta
1.6     ! takayama  142: 行和, 列和のリスト. 成分はすべて正であること.
1.1       takayama  143: @item p
1.6     ! takayama  144: 二元分割表のセルの確率のリスト
1.1       takayama  145: @end table
                    146:
1.6     ! takayama  147: @comment --- ここで関数の詳しい説明 ---
        !           148: @comment --- @itemize〜@end itemize は箇条書き ---
        !           149: @comment --- @bullet は黒点付き ---
1.1       takayama  150: @itemize @bullet
                    151: @item
1.6     ! takayama  152: gmvector は Gauss-Manin vector の略である [GM2016].
1.1       takayama  153: @item
1.6     ! takayama  154: gmvector の戻り値は
        !           155: [GM2016] の 6章 p.23 のベクトル Sである.
        !           156: これは
        !           157: [GM2016] の4章で定義されているベクトル F の定数倍であり,
        !           158: その定数は
        !           159: 第一成分が [GM2016] の6章で定義されている級数 S の値と等しく
        !           160: なるように決められている.
1.1       takayama  161: @item
1.6     ! takayama  162:  r1 x r2 分割表を考える.
        !           163:  m+1=r1, n+1=r2 とおく.
        !           164:  正規化定数 Z は分割表 u を (m+1) × (n+1) 行列とするとき p^u/u! の和である.
        !           165:  ここで和は行和列和が @var{beta} であるような u 全体でとる
1.1       takayama  166:  [TKT2015], [GM2016].
1.6     ! takayama  167:  S はこの多項式 Z の p を
1.1       takayama  168: @verbatim
                    169:   [[1,y11,...,y1n],
                    170:    [1,y21,...,y2n],...,
                    171:    [1,ym1, ...,ymn],
                    172:    [1,1, ..., 1]]
                    173: @end verbatim
1.6     ! takayama  174:  (1 が L 字型に並ぶ),
        !           175: と正規化した級数である.
1.1       takayama  176: @item
1.6     ! takayama  177: 2x(n+1)分割表で, gmvector の戻り値を Lauricella  F_D で書くことが
        !           178: 以下のようにできる
        !           179: (b[2][1]-b[1][1] >= 0 の場合).
        !           180: ここで b[1][1], b[1][2] は, それぞれ 1 行目の行和, 2 行目の行和,
        !           181: b[2][i] は i 列目の列和である.
1.1       takayama  182: @comment ekn/Talks/2015-12-3-goto.tex
                    183: @verbatim
                    184: S=F_D(-b[1,1], [-b[2,2],...,-b[2,n+1]], b[2,1]-b[1,1]+1 ; y)/C,
                    185: @end verbatim
                    186: C=b[1,1]! b[2,2]! ... b[2][n+1]! (b[2,1]-b[1,1])!
1.6     ! takayama  187: とおく.
        !           188: 1/C は L 字型の分割表
1.1       takayama  189: @verbatim
                    190: [[b[1,1],       0,      ..., 0       ],
                    191:  [b[2,1]-b[1,1],b[2,2], ..., b[2,n+1]]]
                    192: @end verbatim
1.6     ! takayama  193: に対応.
        !           194: gmvector は
1.1       takayama  195: @verbatim
                    196: [S,(y11/a2) d_11 S,(y12/a3) d_12 S, ..., (y1n/a_(n+1)) d_1n S]
                    197: @end verbatim
1.6     ! takayama  198: である.
        !           199: ここで d_ij は yij についての微分,
1.1       takayama  200: @verbatim
                    201:   [a0,     a1, ...                      ,a_(n+2)]
                    202: = [-b[1,2],-b[1,1],b[2,2], ..., b[2,n+1],b[2,1]]
                    203: @end verbatim
1.6     ! takayama  204: である.
1.1       takayama  205: @item
1.6     ! takayama  206: 周辺和 @var{beta}の時の正規化定数のセル確率 @var{p} に対する値は 多項式に退化した E(k,n) の値で表現できる. 文献 [TKT2015], [GM2016] 参照.
1.1       takayama  207: @item
1.6     ! takayama  208: option crt=1 (crt = Chinese remainder theorem) を与えると, 分散計算をおこなう
1.1       takayama  209: [T2016].
1.6     ! takayama  210: 分散計算用の各種パラメータの設定は
        !           211: gtt_ekn.setup で行なう.
1.1       takayama  212: @end itemize
                    213:
1.6     ! takayama  214: @comment --- @example〜@end example は実行例の表示 ---
        !           215: 例: 次は2 x 2 分割表で行和が [5,1],  列和が [3,3], 各セルの確率が
        !           216: [[1/2,1/3],[1/7,1/5]] の場合の gmvector の値である.
1.1       takayama  217: @example
                    218: [3000] load("gtt_ekn.rr");
                    219: [3001] ekn_gtt.gmvector([[5,1],[3,3]],[[1/2,1/3],[1/7,1/5]])
                    220: [775/27783]
                    221: [200/9261]
                    222: @end example
                    223:
1.6     ! takayama  224: 参考: 2 x m 分割表(Lauricella FD)についてはパッケージ tk_fd でも下記のように同等な
        !           225: 計算ができる.
        !           226: 守備範囲の異なるプログラム同士の比較, debug 用参考.
1.1       takayama  227: @example
                    228: [3080] import("tk_fd.rr");
                    229: [3081] A=tk_fd.marginal2abc([4,5],[2,4,3]);
1.6     ! takayama  230: [-4,[-4,-3],-1]  // 2変数 FD のパラメータ. a,[b1,b2],c
1.1       takayama  231: [3082] tk_fd.fd_hessian2(A[0],A[1],A[2],[1/2,1/3]);
                    232: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                    233: [4483/124416,[ 1961/15552 185/1728 ],
                    234:  [ 79/288 259/864 ]
                    235:  [ 259/864 47/288 ]]
1.6     ! takayama  236: // 戻値は [F=F_D, gradient(F), Hessian(F)]
1.1       takayama  237:
1.6     ! takayama  238: // ekn_gt での例と同じパラメータ.
1.1       takayama  239: [3543] A=tk_fd.marginal2abc([5,1],[3,3]);
                    240: [-5,[-3],-1]
                    241: [3544] tk_fd.fd_hessian2(A[0],A[1],A[2],[(1/3)*(1/7)/((1/2)*(1/5))]);
                    242: Computing Dmat(ca) for parameters B=[-3],X=[ 10/21 ]
                    243: [775/27783,[ 20/147 ],[ 17/42 ]]
                    244: @end example
                    245:
1.6     ! takayama  246: 参考: 一般の A 分布の正規化定数についての Hessian の計算は実験的 package ot_hessian_ahg.rr
        !           247: で実装のテストがされている. (これはまだ未完成のテスト版なので出力形式等も将来的には変更される.)
1.1       takayama  248: @example
                    249: import("ot_hgm_ahg.rr");
                    250: import("ot_hessian_ahg.rr");
                    251: def  htest4() @{
                    252:   extern C11_A;
                    253:   extern C11_Beta;
                    254:   Hess=newmat(7,7);
                    255:   A =C11_A;
                    256:   Beta0= [b0,b1,b2,b3];
                    257:   BaseIdx=[4,5,6];
                    258:   X=[x0,x1,x2,x3,x4,x5,x6];
                    259:   for (I=0; I<7; I++) for (J=0; J<7; J++) @{
                    260:     Idx = [I,J];
                    261:     H=hessian_simplify(A,Beta0,X,BaseIdx,Idx);
                    262:     Hess[I][J]=H;
                    263:     printf("[I,J]=%a, Hessian_ij=%a\n",Idx,H);
                    264:   @}
                    265:   return(Hess);
                    266: @}
                    267: [2917] C11_A;
                    268: [[0,0,0,1,1,1,1],[1,0,0,1,0,1,0],[0,1,1,0,1,0,1],[1,1,0,1,1,0,0]]
                    269: [2918] C11_Beta;
                    270: [166,36,290,214]
                    271: [2919] Ans=htest4$
                    272: [2920] Ans[0][0];
                    273: [[((b1-b0-1)*x4)/(x0^2),[4]],[((b1-b0-1)*x6)/(x0^2),[6]],
                    274:  [(b1^2+(-2*b0-1)*b1+b0^2+b0)/(x0^2),[]],[(x6)/(x0),[6,0]],[(x4)/(x0),[4,0]]]
                    275: @end example
                    276:
1.6     ! takayama  277: @comment --- 参照(リンク)を書く ---
1.1       takayama  278: @table @t
1.6     ! takayama  279: @item 参照
1.1       takayama  280: @ref{gtt_ekn.setup}
                    281: @ref{gtt_ekn.pfaffian_basis}
                    282: @end table
                    283:
1.6     ! takayama  284: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  285: @noindent
                    286: ChangeLog
                    287: @itemize @bullet
                    288: @item
1.6     ! takayama  289:  この関数は
        !           290: [GM2016] のアルゴリズムおよび
        !           291: [T2016] による modular method を用いた高速化を実装したものである.
1.1       takayama  292: @item
1.6     ! takayama  293:  変更を受けたファイルは
1.1       takayama  294:  OpenXM/src/asir-contrib/packages/src/gtt_ekn.rr 1.1, gtt_ekn/ekn_pfaffian_8.rr
                    295: @end itemize
                    296:
                    297:
                    298: @comment **********************************************************
1.6     ! takayama  299: @node gtt_ekn.nc,,, 超幾何関数E(k,n)
1.1       takayama  300: @subsection @code{gtt_ekn.nc}
1.6     ! takayama  301: @comment --- 索引用キーワード
1.1       takayama  302: @findex gtt_ekn.nc
                    303:
                    304: @table @t
                    305: @item gtt_ekn.nc(@var{beta},@var{p})
1.6     ! takayama  306: :: 周辺和 @var{beta}, セルの確率 @var{p} の二元分割表の条件付き確率の正規化定数 Z
        !           307: およびその微分の値を戻す.
1.1       takayama  308: @end table
                    309:
1.6     ! takayama  310: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  311: @table @var
                    312: @item return
1.6     ! takayama  313: ベクトル [Z,[[d_11 Z, d_12 Z, ...], ..., [d_m1 Z, d_m2 Z, ...., d_mn Z]]]
1.1       takayama  314: @item beta
1.6     ! takayama  315: 行和, 列和のリスト. 成分はすべて正であること.
1.1       takayama  316: @item p
1.6     ! takayama  317: 二元分割表のセルの確率のリスト
1.1       takayama  318: @end table
                    319:
1.6     ! takayama  320: @comment --- ここで関数の詳しい説明 ---
        !           321: @comment --- @itemize〜@end itemize は箇条書き ---
        !           322: @comment --- @bullet は黒点付き ---
1.1       takayama  323: @itemize @bullet
                    324: @item
1.6     ! takayama  325:  r1 x r2 分割表を考える.
        !           326:  m=r1, n=r2 とおく.
        !           327:  正規化定数 Z は分割表 u を m × n 行列とするとき p^u/u! の和である.
        !           328:  ここで和は行和列和が @var{beta} であるような u 全体でとる
1.1       takayama  329:  [TKT2015], [GM2016].
1.6     ! takayama  330:  p^u は p_ij^u_ij の積, u! は u_ij! の積である.
        !           331:  d_ij Z で Z の変数 p_ij についての偏微分を表す.
1.1       takayama  332: @item
1.6     ! takayama  333: nc は gmvector の値を元に, [GM2016] の Prop
        !           334:  7.1 に基づいて Z の値を計算する.
1.1       takayama  335: @item
1.6     ! takayama  336: option crt=1 (crt = Chinese remainder theorem) を与えると, 分散計算をおこなう.
        !           337: 分散計算用の各種パラメータの設定は
        !           338: gtt_ekn.setup で行なう.
1.1       takayama  339: @end itemize
                    340:
1.6     ! takayama  341: @comment --- @example〜@end example は実行例の表示 ---
        !           342: 例: 2x3 分割表での Z とその微分の計算.
1.1       takayama  343: @example
                    344: [2237] gtt_ekn.nc([[4,5],[2,4,3]],[[1,1/2,1/3],[1,1,1]]);
                    345: [4483/124416,[ 353/7776 1961/15552 185/1728 ]
                    346: [ 553/20736 1261/15552 1001/13824 ]]
                    347: @end example
                    348:
1.6     ! takayama  349: 参考: 2 x m 分割表(Lauricella FD)についてはパッケージ tk_fd でも下記のように同等な
        !           350: 計算ができる.
1.1       takayama  351: @example
                    352: [3076] import("tk_fd.rr");
                    353: [3077] A=tk_fd.marginal2abc([4,5],[2,4,3]);
                    354: [-4,[-4,-3],-1]
                    355: [3078] tk_fd.ahmat_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                    356: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                    357: [ 1 1 1 ]
                    358: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                    359: [4483/124416,[[353/7776,1961/15552,185/1728],
                    360:               [553/20736,1261/15552,1001/13824]]]
1.6     ! takayama  361: // 戻値は [Z, [[d_11 Z, d_12 Z, d_13 Z],
        !           362: //             [d_21 Z, d_22 Z, d_23 Z]]] の値.
        !           363: //           ここで d_ij は i,j 成分についての微分を表す.
1.1       takayama  364: @end example
                    365:
1.6     ! takayama  366: @comment --- 参照(リンク)を書く ---
1.1       takayama  367: @table @t
1.6     ! takayama  368: @item 参照
1.1       takayama  369: @ref{gtt_ekn.setup}
                    370: @ref{gtt_ekn.lognc}
                    371: @end table
                    372:
1.6     ! takayama  373: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  374: @noindent
                    375: ChangeLog
                    376: @itemize @bullet
                    377: @item
1.6     ! takayama  378:  変更を受けたファイルは
1.1       takayama  379:  OpenXM/src/asir-contrib/packages/src/gtt_ekn.rr 1.1, gtt_ekn/ekn_eval.rr
                    380: @end itemize
                    381:
                    382:
                    383: @comment **********************************************************
1.6     ! takayama  384: @node gtt_ekn.lognc,,, 超幾何関数E(k,n)
1.1       takayama  385: @subsection @code{gtt_ekn.lognc}
1.6     ! takayama  386: @comment --- 索引用キーワード
1.1       takayama  387: @findex gtt_ekn.lognc
                    388:
                    389: @table @t
                    390: @item gtt_ekn.lognc(@var{beta},@var{p})
1.6     ! takayama  391: :: 周辺和 @var{beta}, セルの確率 @var{p} の二元分割表の条件付き確率の正規化定数 Z
        !           392: の log の近似値およびその微分の近似値を戻す.
1.1       takayama  393: @end table
                    394:
1.6     ! takayama  395: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  396: @table @var
                    397: @item return
1.6     ! takayama  398: ベクトル [log(Z), [[d_11 log(Z), d_12 log(Z), ...], [d_21 log(Z),...], ... ]
1.1       takayama  399: @item beta
1.6     ! takayama  400: 行和, 列和のリスト. 成分はすべて正であること.
1.1       takayama  401: @item p
1.6     ! takayama  402: 二元分割表のセルの確率のリスト
1.1       takayama  403: @end table
                    404:
1.6     ! takayama  405: @comment --- ここで関数の詳しい説明 ---
        !           406: @comment --- @itemize〜@end itemize は箇条書き ---
        !           407: @comment --- @bullet は黒点付き ---
1.1       takayama  408: @itemize @bullet
                    409: @item
1.6     ! takayama  410: 条件付き最尤推定に利用する [TKT2015].
        !           411: @item option crt=1 (crt = Chinese remainder theorem) を与えると, 分散計算をおこなう.
        !           412: 分散計算用の各種パラメータの設定は
        !           413: gtt_ekn.setup で行なう.
1.1       takayama  414: @end itemize
                    415:
1.6     ! takayama  416: @comment --- @example〜@end example は実行例の表示 ---
        !           417: 例: 2 × 3 分割表での例. 第一成分のみ近似値.
1.1       takayama  418: @example
                    419: [2238] gtt_ekn.lognc([[4,5],[2,4,3]],[[1,1/2,1/3],[1,1,1]]);
                    420: [-3.32333832422461674630,[ 5648/4483 15688/4483 13320/4483 ]
                    421: [ 3318/4483 10088/4483 9009/4483 ]]
                    422: @end example
                    423:
1.6     ! takayama  424: 参考: 2 x m 分割表(Lauricella FD)についてはパッケージ tk_fd でも下記のように同等な
        !           425: 計算ができる.
1.1       takayama  426: @example
                    427: [3076] import("tk_fd.rr");
                    428: [3077] A=tk_fd.marginal2abc([4,5],[2,4,3]);
                    429: [-4,[-4,-3],-1]
                    430: [3078] tk_fd.log_ahmat_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                    431: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                    432: [ 1 1 1 ]
                    433: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                    434: [-3.32333832422461674639485797719209322217260539267246045320,
                    435:  [[1.2598706, 3.499442, 2.971224],
                    436:   [0.7401293, 2.250278, 2.009591]]]
1.6     ! takayama  437: // 戻値は [log(Z),
1.1       takayama  438: //          [[d_11 log(Z), d_12 log(Z), d_13 log(Z)],
                    439: //           [d_21 log(Z), d_22 log(Z), d_23 log(Z)]]]
1.6     ! takayama  440: // の近似値.
1.1       takayama  441: @end example
                    442:
1.6     ! takayama  443: @comment --- 参照(リンク)を書く ---
1.1       takayama  444: @table @t
1.6     ! takayama  445: @item 参照
1.1       takayama  446: @ref{gtt_ekn.setup}
                    447: @ref{gtt_ekn.nc}
                    448: @end table
                    449:
1.6     ! takayama  450: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  451: @noindent
                    452: ChangeLog
                    453: @itemize @bullet
                    454: @item
1.6     ! takayama  455:  変更を受けたファイルは
1.1       takayama  456:  OpenXM/src/asir-contrib/packages/src/gtt_ekn.rr 1.1.
                    457: @end itemize
                    458:
                    459: @comment **********************************************************
1.6     ! takayama  460: @node gtt_ekn.expectation,,, 超幾何関数E(k,n)
1.1       takayama  461: @subsection @code{gtt_ekn.expectation}
1.6     ! takayama  462: @comment --- 索引用キーワード
1.1       takayama  463: @findex gtt_ekn.expectation
                    464:
                    465: @table @t
                    466: @item gtt_ekn.expectation(@var{beta},@var{p})
1.6     ! takayama  467: :: 周辺和 @var{beta}, セルの確率 @var{p} の二元分割表の期待値を計算する.
1.1       takayama  468: @end table
                    469:
1.6     ! takayama  470: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  471: @table @var
                    472: @item return
1.6     ! takayama  473: 二元分割表の各セルの期待値のリスト.
1.1       takayama  474: @item beta
1.6     ! takayama  475: 行和, 列和のリスト. 成分はすべて正であること.
1.1       takayama  476: @item p
1.6     ! takayama  477: 二元分割表のセルの確率のリスト
1.1       takayama  478: @end table
                    479:
1.6     ! takayama  480: @comment --- ここで関数の詳しい説明 ---
        !           481: @comment --- @itemize〜@end itemize は箇条書き ---
        !           482: @comment --- @bullet は黒点付き ---
1.1       takayama  483: @itemize @bullet
                    484: @item
1.6     ! takayama  485: [GM2016] の Algorithm 7.8 の実装.
        !           486: @item option crt=1 (crt = Chinese remainder theorem) を与えると, 分散計算をおこなう.
        !           487: 分散計算用の各種パラメータの設定は
        !           488: gtt_ekn.setup で行なう.
        !           489: @item option index を与えると, 指定された成分の期待値のみ計算する.
        !           490: たとえば 2 x 2 分割表で index=[[0,0],[1,1]] と指定すると, 1 のある成分の期待値のみ計算する.
1.1       takayama  491: @end itemize
                    492:
1.6     ! takayama  493: @comment --- @example〜@end example は実行例の表示 ---
1.1       takayama  494:
1.6     ! takayama  495: 2×2, 3×3 の分割表の期待値計算例.
1.1       takayama  496: @example
                    497: [2235] gtt_ekn.expectation([[1,4],[2,3]],[[1,1/3],[1,1]]);
                    498: [ 2/3 1/3 ]
                    499: [ 4/3 8/3 ]
                    500: [2236] gtt_ekn.expectation([[4,5],[2,4,3]],[[1,1/2,1/3],[1,1,1]]);
                    501: [ 5648/4483 7844/4483 4440/4483 ]
                    502: [ 3318/4483 10088/4483 9009/4483 ]
                    503:
                    504: [2442] gtt_ekn.expectation([[4,14,9],[11,6,10]],[[1,1/2,1/3],[1,1/5,1/7],[1,1,1]]);
                    505: [ 207017568232262040/147000422096729819 163140751505489940/147000422096729819
                    506:                                         217843368649167296/147000422096729819 ]
                    507: [ 1185482401011137878/147000422096729819 358095302885438604/147000422096729819
                    508:                                          514428205457640984/147000422096729819 ]
                    509: [ 224504673820628091/147000422096729819 360766478189450370/147000422096729819
                    510:                                         737732646860489910/147000422096729819 ]
                    511: @end example
                    512:
1.6     ! takayama  513: 参考: 2 x m 分割表(Lauricella FD)についてはパッケージ tk_fd でも下記のように同等な
        !           514: 計算ができる.
1.1       takayama  515: @example
                    516: [3076] import("tk_fd.rr");
                    517: [3077] A=tk_fd.marginal2abc([4,5],[2,4,3]);
                    518: [-4,[-4,-3],-1]
                    519: [3078] tk_fd.expectation_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                    520: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                    521: [ 1 1 1 ]
                    522: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                    523: [[5648/4483,7844/4483,4440/4483],
                    524:  [3318/4483,10088/4483,9009/4483]]
1.6     ! takayama  525: // 各セルの期待値.
1.1       takayama  526: @end example
                    527:
1.6     ! takayama  528: 参考: 一般の A 分布の計算は ot_hgm_ahg.rr. まだ実験的なため, module 化されていない.
        !           529: ot_hgm_ahg.rr についての参考文献:
1.1       takayama  530: K.Ohara, N.Takayama, Pfaffian Systems of A-Hypergeometric Systems II --- Holonomic Gradient Method, arxiv:1505.02947
                    531: @example
                    532: [3237] import("ot_hgm_ahg.rr");
1.6     ! takayama  533: // 2 x 2 分割表.
1.1       takayama  534: [3238] hgm_ahg_expected_values_contiguity([[0,0,1,1],[1,0,1,0],[0,1,0,1]],
                    535:         [9,6,8],[1/2,1/3,1/5,1/7],[x1,x2,x3,x4]|geometric=1);
                    536: oohg_native=0, oohg_curl=1
                    537: [1376777025/625400597,1750225960/625400597,
                    538:  2375626557/625400597,3252978816/625400597]
1.6     ! takayama  539: // 2 x 2 分割表の期待値.
1.1       takayama  540:
1.6     ! takayama  541: // 2 x 3 分割表.
1.1       takayama  542: [3238] hgm_ahg_expected_values_contiguity(
                    543:  [[0,0,0,1,1,1],[1,0,0,1,0,0],[0,1,0,0,1,0],[0,0,1,0,0,1]],
                    544:  [5,2,4,3],[1,1/2,1/3,1,1,1],[x1,x2,x3,x4,x5,x6]|geometric=1);
                    545: [5648/4483,7844/4483,4440/4483,3318/4483,10088/4483,9009/4483]
1.6     ! takayama  546: // 2 x 3 分割表の期待値. 上と同じ問題.
1.1       takayama  547: @end example
                    548:
1.6     ! takayama  549: 3 x 3 分割表. 構造的0が一つ.
1.1       takayama  550: @example
                    551: /*
1.6     ! takayama  552:   dojo, p.221 のデータ.  成績3以下の生徒は集めてひとつに.
1.1       takayama  553:   2 1 1
                    554:   8 3 3
                    555:   0 2 6
                    556:
                    557:   row sum: 4,14,8
                    558:   column sum: 10,6,10
1.6     ! takayama  559:   0 を一つ含むので, (3,6) 型の A から 7 列目を抜く.
1.1       takayama  560: */
                    561:
                    562: A=[[0,0,0,1,1,1, 0,0],
                    563:    [0,0,0,0,0,0, 1,1],
                    564:    [1,0,0,1,0,0, 0,0],
                    565:    [0,1,0,0,1,0, 1,0],
                    566:    [0,0,1,0,0,1, 0,1]];
                    567: B=[14,8,10,6,10];
                    568: hgm_ahg_expected_values_contiguity(A,B,[1,1/2,1/3,1,1/5,1/7,1,1],
1.6     ! takayama  569:                 [x1,x2,x3,x4,x5,x6,x7,x8]|geometric=1);
1.1       takayama  570:
1.6     ! takayama  571: // 答.
1.1       takayama  572: [14449864949304/9556267369631,
                    573:  10262588586540/9556267369631, 13512615942680/9556267369631,
                    574:  81112808747006/9556267369631,
                    575:  21816297744346/9556267369631, 30858636683482/9556267369631,
                    576:
                    577:  25258717886900/9556267369631,51191421070148/9556267369631]
                    578: @end example
                    579:
1.6     ! takayama  580: 3 x 3 分割表.
1.1       takayama  581: @example
                    582: /*
1.6     ! takayama  583:  上のデータで 0 を 1 に変更.
1.1       takayama  584:   2 1 1
                    585:   8 3 3
                    586:   1 2 6
                    587:
                    588:   row sum: 4,14,9
                    589:   column sum: 11,6,10
                    590: */
                    591: A=[[0,0,0,1,1,1,0,0,0],
                    592:    [0,0,0,0,0,0,1,1,1],
                    593:    [1,0,0,1,0,0,1,0,0],
                    594:    [0,1,0,0,1,0,0,1,0],
                    595:    [0,0,1,0,0,1,0,0,1]];
                    596: B=[14,9,11,6,10];
                    597: hgm_ahg_expected_values_contiguity(A,B,[1,1/2,1/3,1,1/5,1/7,1,1,1],
                    598:                               [x1,x2,x3,x4,x5,x6,x7,x8]|geometric=1);
                    599:
1.6     ! takayama  600: // 期待値, 答.   x9 を指定していないので, 9番目の期待値は出力してない.
1.1       takayama  601: [207017568232262040/147000422096729819,
                    602:  163140751505489940/147000422096729819,217843368649167296/147000422096729819,
                    603:  1185482401011137878/147000422096729819,
                    604:  358095302885438604/147000422096729819,514428205457640984/147000422096729819,
                    605:  224504673820628091/147000422096729819,360766478189450370/147000422096729819]
                    606:
1.6     ! takayama  607: // Z やその微分の計算は hgm_ahg_contiguity 関数がおこなうが, これの簡易インターフェースは
        !           608: // まだ書いてない.
1.1       takayama  609: @end example
                    610:
                    611:
                    612:
1.6     ! takayama  613: @comment --- 参照(リンク)を書く ---
1.1       takayama  614: @table @t
1.6     ! takayama  615: @item 参照
1.1       takayama  616: @ref{gtt_ekn.setup}
                    617: @ref{gtt_ekn.nc}
                    618: @end table
                    619:
1.6     ! takayama  620: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  621: @noindent
                    622: ChangeLog
                    623: @itemize @bullet
                    624: @item
1.6     ! takayama  625:  変更を受けたファイルは
1.1       takayama  626:  OpenXM/src/asir-contrib/packages/src/gtt_ekn.rr 1.1.
                    627: @end itemize
                    628:
                    629:
                    630: @comment **********************************************************
1.6     ! takayama  631: @comment --- ◯◯◯◯  の説明
        !           632: @comment --- 個々の関数の説明の開始 ---
        !           633: @comment --- section 名を正確に ---
        !           634: @node gtt_ekn.setup,,, 超幾何関数E(k,n)
1.1       takayama  635: @subsection @code{gtt_ekn.setup}
1.6     ! takayama  636: @comment --- 索引用キーワード
1.1       takayama  637: @findex gtt_ekn.setup
                    638:
                    639: @table @t
                    640: @item gtt_ekn.setup()
1.6     ! takayama  641: :: 分散計算用の環境設定をおこなう. 現在の環境を報告する.
1.1       takayama  642: @end table
                    643:
1.6     ! takayama  644: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  645: @table @var
                    646: @item return
                    647:
                    648: @end table
                    649:
1.6     ! takayama  650: @comment --- ここで関数の詳しい説明 ---
        !           651: @comment --- @itemize〜@end itemize は箇条書き ---
        !           652: @comment --- @bullet は黒点付き ---
1.3       takayama  653: @itemize @bullet
1.6     ! takayama  654: @item 使用するプロセスと素数の個数, 最小の素数を表示する. 準備されていない場合はその旨を表示.
        !           655: @item このパッケージでの分散計算は複数のcpuを搭載した計算機で実行されることを想定している.
        !           656: @item option nps (または number_of_processes)を与えると指定した数だけプロセスを用意する.
        !           657: @item option nprm (または number_of_primes)を与えるとnprmが文字列の場合指定された素数リストのファイルを読み込む. nprmが自然数の場合さらにoption minp (minp =MINimum Prime)を与えるとminpより大きな素数をnprm個生成する. その際option fgp (または file_of_generated_primes)を与えると生成した素数リストをファイル名をfgpとして保存する.
        !           658: @item 上記のoption を指定しなかった場合次のデフォルト値が用いられる. nps=1. nprm=10. fgp=0.
        !           659: @item option report=1を与えると現在の環境の報告のみを行う. setup(|report=1)の別名としてreport関数を使用することもできる.
        !           660: @item option subprogs=[file1,file2,...] により分散計算の子供プロセスにロードすべきファイル file1, file2, ... を指定する. default は subprogs=["gtt_ekn/childprocess.rr"] である.
        !           661: @item gtt_ekn.set_debug_mode(Mode) で Ekn_debug の値を設定する.
1.1       takayama  662: @end itemize
                    663:
1.6     ! takayama  664: @comment --- @example〜@end example は実行例の表示 ---
        !           665: 例: 素数のリストを生成してファイル p.txt へ書き出す.
1.1       takayama  666: @example
                    667: gtt_ekn.setup(|nps=2,nprm=20,minp=10^10,fgp="p.txt")$
                    668: @end example
                    669:
                    670:
1.6     ! takayama  671: @comment --- 参照(リンク)を書く ---
1.1       takayama  672: @table @t
1.6     ! takayama  673: @item 参照
1.1       takayama  674: @ref{gtt_ekn.nc}
                    675: @ref{gtt_ekn.gmvector}
                    676: @end table
                    677:
1.6     ! takayama  678: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  679: @noindent
                    680: ChangeLog
                    681: @itemize @bullet
                    682: @item
1.6     ! takayama  683:  変更を受けたファイルは
1.1       takayama  684:  OpenXM/src/asir-contrib/packages/src/gtt_ekn.rr 1.1,
                    685:  gtt_ekn/g_mat_fac.rr
                    686:
                    687: @end itemize
                    688:
                    689: @comment **********************************************************
1.6     ! takayama  690: @comment --- ◯◯◯◯  の説明
        !           691: @comment --- 個々の関数の説明の開始 ---
        !           692: @comment --- section 名を正確に ---
        !           693: @node gtt_ekn.upAlpha,,, 超幾何関数E(k,n)
1.1       takayama  694: @subsection @code{gtt_ekn.upAlpha}
1.6     ! takayama  695: @comment --- 索引用キーワード
1.1       takayama  696: @findex gtt_ekn.upAlpha
                    697:
                    698: @table @t
                    699: @item gtt_ekn.upAlpha(@var{i},@var{k},@var{n})
                    700: ::
                    701: @end table
                    702:
1.6     ! takayama  703: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.1       takayama  704: @table @var
1.6     ! takayama  705: @item i  a_i を a_i+1 と変化させる contiguity relation.
        !           706: @item k  E(k+1,n+k+2)型の超幾何関数の k. 分割表では (k+1)×(n+1).
        !           707: @item n  E(k+1,n+k+2)型の超幾何関数の n. 分割表では (k+1)×(n+1).
        !           708: @item return  contiguity relation の pfaffian_basis についての行列表現を戻す. [GM2016] の Cor 6.3.
1.1       takayama  709: @end table
                    710:
1.6     ! takayama  711: @comment --- ここで関数の詳しい説明 ---
        !           712: @comment --- @itemize〜@end itemize は箇条書き ---
        !           713: @comment --- @bullet は黒点付き ---
1.1       takayama  714: @itemize @bullet
                    715: @item
1.6     ! takayama  716:  upAlpha は [GM2016] の Cor 6.3 の行列 U_i を戻す.
        !           717: @item 関連する各関数の簡潔な説明と例も加える.
        !           718: @item a_i を a_i-1 と変化させたい場合は関数 downAlpha を用いる.
        !           719: @item a_i と分割表の周辺和を見るには, 関数 marginaltoAlpha([行和,列和]) を用いる.
1.1       takayama  720: @item
1.6     ! takayama  721:    pfaffian_basis は [GM2016] の4章のベクトル F に対応する偏微分を戻す.
1.1       takayama  722: @end itemize
                    723:
1.6     ! takayama  724: @comment --- @example〜@end example は実行例の表示 ---
        !           725: 例: 以下の例は 2×2分割表(E(2,4)), 2×3分割表(E(2,5))の場合である.
        !           726: [2225] までは出力を略している.
1.1       takayama  727: @example
                    728: [2221] gtt_ekn.marginaltoAlpha([[1,4],[2,3]]);
                    729: [[a_0,-4],[a_1,-1],[a_2,3],[a_3,2]]
1.6     ! takayama  730: [2222] gtt_ekn.upAlpha(1,1,1);  // E(2,4) の a_1 方向の
        !           731:                                 //     contiguity を表現する行列
        !           732: [2223] gtt_ekn.upAlpha(2,1,1);  // E(2,4) の a_2 方向
        !           733: [2224] gtt_ekn.upAlpha(3,1,1);  // E(2,4) の a_3 方向
1.1       takayama  734: [2225] function f(x_1_1);
                    735: [2232] gtt_ekn.pfaffian_basis(f(x_1_1),1,1);
                    736: [ f(x_1_1) ]
                    737: [ (f{1}(x_1_1)*x_1_1)/(a_2) ]
                    738: [2233] function f(x_1_1,x_1_2);
                    739: f() redefined.
1.6     ! takayama  740: [2234] gtt_ekn.pfaffian_basis(f(x_1_1,x_1_2),1,2); // E(2,5), 2*3 分割表
1.1       takayama  741: [ f(x_1_1,x_1_2) ]
                    742: [ (f{1,0}(x_1_1,x_1_2)*x_1_1)/(a_2) ]
                    743: [ (f{0,1}(x_1_1,x_1_2)*x_1_2)/(a_3) ]
                    744: @end example
                    745:
                    746:
1.6     ! takayama  747: @comment --- 参照(リンク)を書く ---
1.1       takayama  748: @table @t
1.6     ! takayama  749: @item 参照
1.1       takayama  750: @ref{gtt_ekn.nc}
                    751: @ref{gtt_ekn.gmvector}
                    752: @end table
                    753:
1.6     ! takayama  754: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.1       takayama  755: @noindent
                    756: ChangeLog
                    757: @itemize @bullet
                    758: @item
1.6     ! takayama  759:  この関数は [GM2016]
        !           760: で与えられたアルゴリズムに従い contiguity relation を導出する.
1.1       takayama  761: @item
1.6     ! takayama  762:  変更を受けたファイルは
1.1       takayama  763:  OpenXM/src/asir-contrib/packages/src/gtt_ekn/ekn_pfaffian_8.rr 1.1.
                    764: @end itemize
                    765:
                    766:
1.5       takayama  767: @comment **********************************************************
1.6     ! takayama  768: @comment --- ◯◯◯◯  の説明
        !           769: @comment --- 個々の関数の説明の開始 ---
        !           770: @comment --- section 名を正確に ---
        !           771: @node gtt_ekn.cmle,,, 超幾何関数E(k,n)
1.5       takayama  772: @subsection @code{gtt_ekn.cmle}
1.6     ! takayama  773: @comment --- 索引用キーワード
1.5       takayama  774: @findex gtt_ekn.cmle
                    775:
                    776: @table @t
1.6     ! takayama  777: @item gtt_ekn.cmle(@var{u}) u を観測データとするとき, P(U=u | row sum, column sum = these of U) を最大化する, 各セルの確率の近似値を求める.
1.5       takayama  778: ::
                    779: @end table
                    780:
1.6     ! takayama  781: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.5       takayama  782: @table @var
1.6     ! takayama  783: @item u  観測データ(分割表)
        !           784: @item return  セルの確率(分割表形式)
1.5       takayama  785: @end table
                    786:
1.6     ! takayama  787: @comment --- ここで関数の詳しい説明 ---
        !           788: @comment --- @itemize〜@end itemize は箇条書き ---
        !           789: @comment --- @bullet は黒点付き ---
1.5       takayama  790: @itemize @bullet
1.6     ! takayama  791: @item u を観測データとするとき, P(U=u | row sum, column sum = these of U) を最大化する, 各セルの確率の近似値を求める.
        !           792: @item optional parameter で algorithm の振る舞い(たとえば有理数を近似して, 分母分子が小さい有理数にする, gradient descent の step幅)を調整すべきだが, これは作業中. 2017.03.03
1.5       takayama  793: @end itemize
                    794:
1.6     ! takayama  795: @comment --- @example〜@end example は実行例の表示 ---
        !           796: 例: 2 x 4 分割表.
1.5       takayama  797: @example
                    798: U=[[1,1,2,3],[1,3,1,1]];
                    799: gtt_ekn.cmle(U);
                    800:  [[ 1 1 2 3 ]
                    801:   [ 1 3 1 1 ],[[7,6],[2,4,3,4]],   // Data, row sum, column sum
                    802:  [ 1 67147/183792 120403/64148 48801/17869 ]  // probability obtained.
                    803:  [ 1 1 1 1 ]]
                    804: @end example
                    805:
1.6     ! takayama  806: 例: 上の例は次の関数に.
1.5       takayama  807: @example
                    808: gtt_ekn.cmle_test3();
                    809: @end example
                    810:
1.6     ! takayama  811: @comment --- 参照(リンク)を書く ---
1.5       takayama  812: @table @t
1.6     ! takayama  813: @item 参照
1.5       takayama  814: @ref{gtt_ekn.expectation}
                    815: @end table
                    816:
1.6     ! takayama  817: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.5       takayama  818: @noindent
                    819: ChangeLog
                    820: @itemize @bullet
1.6     ! takayama  821: @item  gtt_ekn/mle.rr に本体がある.
        !           822: @item  gtt_ekn.rr の cmle 関数は wrapper.
1.5       takayama  823: @end itemize
                    824: @comment end cmle.
                    825:
                    826:
                    827:
1.6     ! takayama  828: @node modular計算,,, 2元分割表HGMの関数
        !           829: @chapter modular計算
1.4       takayama  830:
                    831: @menu
                    832: * gtt_ekn.chinese_itor::
                    833: @end menu
                    834:
1.6     ! takayama  835: @node 中国剰余定理とitor,,, modular計算
        !           836: @section 中国剰余定理とitor
1.4       takayama  837:
                    838: @comment **********************************************************
1.6     ! takayama  839: @comment --- ◯◯◯◯  の説明
        !           840: @comment --- 個々の関数の説明の開始 ---
        !           841: @comment --- section 名を正確に ---
1.4       takayama  842: @node gtt_ekn.chinese_itor,,,
                    843: @subsection @code{gtt_ekn.chinese_itor}
1.6     ! takayama  844: @comment --- 索引用キーワード
        !           845: @findex gtt_ekn.chinese_itor 中国剰余定理とitor
1.4       takayama  846:
                    847: @table @t
                    848: @item gtt_ekn.chinese_itor(@var{data},@var{idlist})
1.6     ! takayama  849: :: mod p で計算した結果(ベクトル)から chinese remainder theorem, itor(integer to rational) で有理数ベクトルを得る.
1.4       takayama  850: @end table
                    851:
1.6     ! takayama  852: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.4       takayama  853: @table @var
1.6     ! takayama  854: @item return  [val, n]  ここで val は答え. また, n = n1*n2*...
        !           855: @item data   [[val1,n1],[val2,n2], ...], ここで val mod n1 = val1, val mod n2 = val2,...
        !           856: @item idlist  chinese, itor を実行するサーバIDのリスト.
1.4       takayama  857: @end table
                    858:
1.6     ! takayama  859: @comment --- ここで関数の詳しい説明 ---
        !           860: @comment --- @itemize〜@end itemize は箇条書き ---
        !           861: @comment --- @bullet は黒点付き ---
1.4       takayama  862: @itemize @bullet
1.6     ! takayama  863: @item 中国剰余定理を用いて val0 mod n1 = val1, val0 mod n2 = val2, ... となる val0 を求める. val に algorithm itor を適用する.
        !           864: @item sqrt(n) より val0 が大きい時は itor が適用されて val0 が有理数 val=a/b に変換される. つまり b*x =1 mod n となる逆数 x を考えて, x*a % n = val0 となる数 val を戻す. 見つからないときは failure を戻す.
1.4       takayama  865: @end itemize
                    866:
1.6     ! takayama  867: @comment --- @example〜@end example は実行例の表示 ---
        !           868: 例: [3!, 5^3*3!]=[6,750] が戻り値.
        !           869: 6 mod 109 =6, 750 mod 109=96 が最初の引数の [[6,96],109]. 以下同様.
1.4       takayama  870: @example
                    871: gtt_ekn.setup(|nps=2,nprm=3,minp=101,fgp="p_small.txt");
                    872: SS=gtt_ekn.get_svalue();
                    873: SS[0];
                    874:   [103,107,109]   // list of primes
                    875: SS[1];
                    876:   [0,2]           // list of server ID's
                    877: gtt_ekn.chinese_itor([[[ 6,96 ],109],[[ 6,29 ],103],[[ 6,1 ],107]],SS[1]);
                    878:   [[ 6 750 ],1201289]
                    879:
1.6     ! takayama  880: // 引数はスカラーでもよい.
1.4       takayama  881: gtt_ekn.chinese_itor([[96,109],[29,103]],SS[1]);
                    882:   [[ 750 ],11227]
                    883: @end example
                    884:
                    885:
1.6     ! takayama  886: @comment --- @example〜@end example は実行例の表示 ---
        !           887: 例: gtt_ekn/childprocess.rr (server で実行される) の関数 chinese (chinese remainder theorem) と euclid.
1.4       takayama  888: @example
                    889: load("gtt_ekn/childprocess.rr");
                    890: chinese([newvect(2,[6,29]),103],[newvect(2,[6,750]),107*109]);
1.6     ! takayama  891:   // mod 103 で [6,29], mod (107*109) で [6,750] となる数を mod 103*(107*109)
        !           892:   // で求めると,
1.4       takayama  893:   [[ 6 750 ],1201289]
1.6     ! takayama  894: euclid(3,103);  // mod 103 での 3 の逆数. つまり 1/3
1.4       takayama  895:   -34
1.6     ! takayama  896: 3*(-34) % 103; // 確かに逆数.
1.4       takayama  897:    1
                    898: @end example
                    899:
1.6     ! takayama  900: @comment --- @example〜@end example は実行例の表示 ---
        !           901: 例: gtt_ekn/childprocess.rr (server で実行される) の関数 itor (integer to rational) の例.
        !           902: itor(Y,Q,Q2,Idx) では Y < Q2 なら Y がそのまま戻る.  Idx は 内部用の index で好きな数でよい. 戻り値の第2成分となる.
1.4       takayama  903: @example
                    904: load("gtt_ekn/childprocess.rr");
                    905: for (I=1;I<11; I++) print([I,itor(I,11,3,0)]);
                    906: [1,[1,0]]
                    907: [2,[2,0]]
1.6     ! takayama  908: [3,[-2/3,0]] //euclid(3,11); ->4,  4*(-2)%11 -> 3 なので確かに -2/3 は元の数の候補
1.4       takayama  909: [4,[failure,0]]
                    910: [5,[-1/2,0]]
                    911: [6,[1/2,0]]
                    912: [7,[-1/3,0]]
                    913: [8,[failure,0]]
                    914: [9,[-2,0]]
                    915: [10,[-1,0]]
                    916: @end example
                    917:
                    918:
1.6     ! takayama  919: @comment --- 参照(リンク)を書く ---
1.4       takayama  920: @table @t
1.6     ! takayama  921: @item 参照
1.4       takayama  922: @ref{gtt_ekn.setup}
                    923: @end table
                    924:
1.6     ! takayama  925: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.4       takayama  926: @noindent
                    927: ChangeLog
                    928: @itemize @bullet
                    929: @item
1.6     ! takayama  930:  関連ファイルは
1.4       takayama  931:  gtt_ekn/g_mat_fac.rr
                    932:  gtt_ekn/childprocess.rr
                    933: @end itemize
                    934:
                    935:
1.1       takayama  936:
1.6     ! takayama  937: @comment --- おまじない ---
1.1       takayama  938: @node Index,,, Top
                    939: @unnumbered Index
                    940: @printindex fn
                    941: @printindex cp
                    942: @iftex
                    943: @vfill @eject
                    944: @end iftex
                    945: @summarycontents
                    946: @contents
                    947: @bye
1.6     ! takayama  948: @comment --- おまじない終り ---
1.1       takayama  949:
                    950:
1.6     ! takayama  951: @comment テンプレート.  start_of_template.
1.5       takayama  952: @comment **********************************************************
1.6     ! takayama  953: @comment --- ◯◯◯◯  の説明
        !           954: @comment --- 個々の関数の説明の開始 ---
        !           955: @comment --- section 名を正確に ---
        !           956: @node gtt_ekn.hoge,,, 超幾何関数E(k,n)
1.5       takayama  957: @subsection @code{gtt_ekn.hoge}
1.6     ! takayama  958: @comment --- 索引用キーワード
1.5       takayama  959: @findex gtt_ekn.hoge
                    960:
                    961: @table @t
                    962: @item gtt_ekn.hoge(@var{i})
                    963: ::
                    964: @end table
                    965:
1.6     ! takayama  966: @comment --- 引数の簡単な説明 ---  以下まだ書いてない.
1.5       takayama  967: @table @var
                    968: @item i  hage
                    969: @item return
                    970: @end table
                    971:
1.6     ! takayama  972: @comment --- ここで関数の詳しい説明 ---
        !           973: @comment --- @itemize〜@end itemize は箇条書き ---
        !           974: @comment --- @bullet は黒点付き ---
1.5       takayama  975: @itemize @bullet
1.6     ! takayama  976: @item 説明.
1.5       takayama  977: @end itemize
                    978:
1.6     ! takayama  979: @comment --- @example〜@end example は実行例の表示 ---
        !           980: 例:
1.5       takayama  981: @example
                    982: [2221] gtt_ekn.hoge([[1,4],[2,3]]);
                    983: @end example
                    984:
                    985:
1.6     ! takayama  986: @comment --- 参照(リンク)を書く ---
1.5       takayama  987: @table @t
1.6     ! takayama  988: @item 参照
1.5       takayama  989: @ref{gtt_ekn.nc}
                    990: @ref{gtt_ekn.gmvector}
                    991: @end table
                    992:
1.6     ! takayama  993: @comment --- ChangeLog を書く. ソースコードの位置. 変更日時 など CVSサーバを見るため
1.5       takayama  994: @noindent
                    995: ChangeLog
                    996: @itemize @bullet
                    997: @item
                    998: @end itemize
                    999: @comment end_of_template
                   1000:
                   1001:
1.6     ! takayama 1002: // 2 x m 分割表において似た機能を有する関数の利用例を参考までに記載する;
        !          1003: // 正規化定数とその微分関連.
        !          1004: // その1.
1.1       takayama 1005: [3077] A=tk_fd.marginal2abc([4,5],[2,4,3]);
                   1006: [-4,[-4,-3],-1]
                   1007: [3078] tk_fd.ahmat_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                   1008: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                   1009: [ 1 1 1 ]
                   1010: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                   1011: [4483/124416,[[353/7776,1961/15552,185/1728],[553/20736,1261/15552,1001/13824]]]
1.6     ! takayama 1012: // 戻値は [Z, [[d_11 Z, d_12 Z, d_13 Z],[d_21 Z, d_22 Z, d_23 Z]]] の値.
1.1       takayama 1013:
1.6     ! takayama 1014: // その2.
1.1       takayama 1015: [3079] tk_fd.log_ahmat_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                   1016: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                   1017: [ 1 1 1 ]
                   1018: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                   1019: [-3.32333832422461674639485797719209322217260539267246045320,
                   1020:  [[1.25987062235110417131385233102832924380994869507026544724,3.49944233772027660049074280615659156814633058219942003122,2.97122462636627258532232879768012491635065804149007361142],
                   1021:   [0.740129377648895828686147668971670756190051304929734552754,2.25027883113986169975462859692170421592683470890028998438,2.00959179121124247155922373410662502788311398616997546285]]]
1.6     ! takayama 1022: // 戻値は [log(Z),
1.1       takayama 1023: //          [[d_11 log(Z), d_12 log(Z), d_13 log(Z)],
                   1024: //           [d_21 log(Z), d_22 log(Z), d_23 log(Z)]]]
1.6     ! takayama 1025: // の近似値.
1.1       takayama 1026:
1.6     ! takayama 1027: // その3.
1.1       takayama 1028: [3082] fd_hessian2(A[0],A[1],A[2],[1/2,1/3]);
                   1029: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                   1030: [4483/124416,[ 1961/15552 185/1728 ],
                   1031:  [ 79/288 259/864 ]
                   1032:  [ 259/864 47/288 ]]
1.6     ! takayama 1033: // 戻値は [F=F_D, gradient(F), Hessian(F)]
1.1       takayama 1034:
1.6     ! takayama 1035: // 参考.
        !          1036: // ygahvec で巾関数分の調整. 独立した関数はないようだ.
1.1       takayama 1037:
                   1038: //-----------------------------------------------------------------------
1.6     ! takayama 1039: // 2 x m 分割表において似た機能を有する関数の利用例を参考までに記載する;
        !          1040: // 期待値関連.
1.1       takayama 1041: [3079] A=tk_fd.marginal2abc([4,5],[2,4,3]);
                   1042: [-4,[-4,-3],-1]
                   1043: [3080] tk_fd.expectation_abc(A[0],A[1],A[2],[[1,1/2,1/3],[1,1,1]]);
                   1044: RS=[ 4 5 ], CSnew=[ 2 4 3 ], Ynew=[ 1 1/2 1/3 ]
                   1045: [ 1 1 1 ]
                   1046: Computing Dmat(ca) for parameters B=[-4,-3],X=[ 1/2 1/3 ]
                   1047: [[5648/4483,7844/4483,4440/4483],
                   1048:  [3318/4483,10088/4483,9009/4483]]
1.6     ! takayama 1049: // 各セルの期待値.
1.1       takayama 1050:
                   1051: //-----------------------------------------------------------------------
1.6     ! takayama 1052: // ot_hgm_ahg.rr の例.  実験的なため module 化されていない.
1.1       takayama 1053: [3237] import("ot_hgm_ahg.rr");
1.6     ! takayama 1054: // 2 x 2 分割表.
1.1       takayama 1055: [3238] hgm_ahg_expected_values_contiguity([[0,0,1,1],[1,0,1,0],[0,1,0,1]],
                   1056:         [9,6,8],[1/2,1/3,1/5,1/7],[x1,x2,x3,x4]|geometric=1);
                   1057: oohg_native=0, oohg_curl=1
                   1058: [1376777025/625400597,1750225960/625400597,2375626557/625400597,3252978816/625400597]
1.6     ! takayama 1059: // 2 x 2 分割表の期待値.
1.1       takayama 1060:
1.6     ! takayama 1061: // 2 x 3 分割表.
1.1       takayama 1062: [3238] hgm_ahg_expected_values_contiguity(
                   1063:  [[0,0,0,1,1,1],[1,0,0,1,0,0],[0,1,0,0,1,0],[0,0,1,0,0,1]],
                   1064:  [5,2,4,3],[1,1/2,1/3,1,1,1],[x1,x2,x3,x4,x5,x6]|geometric=1);
                   1065: [5648/4483,7844/4483,4440/4483,3318/4483,10088/4483,9009/4483]
1.6     ! takayama 1066: // 2 x 3 分割表の期待値. 上と同じ問題.
1.1       takayama 1067:
                   1068: /*
1.6     ! takayama 1069:   dojo, p.221.  成績3以下の生徒は集めてひとつに.
1.1       takayama 1070:   2 1 1
                   1071:   8 3 3
                   1072:   0 2 6
                   1073:
                   1074:   row sum: 4,14,8
                   1075:   column sum: 10,6,10
1.6     ! takayama 1076:   0 を一つ含むので, (3,6) 型の A から 7 列目を抜く.
1.1       takayama 1077: */
1.6     ! takayama 1078: // 3 x 3 分割表. 構造的0が一つ.
1.1       takayama 1079:
                   1080: A=[[0,0,0,1,1,1, 0,0],
                   1081:    [0,0,0,0,0,0, 1,1],
                   1082:    [1,0,0,1,0,0, 0,0],
                   1083:    [0,1,0,0,1,0, 1,0],
                   1084:    [0,0,1,0,0,1, 0,1]];
                   1085: B=[14,8,10,6,10];
                   1086: hgm_ahg_expected_values_contiguity(A,B,[1,1/2,1/3,1,1/5,1/7,1,1],[x1,x2,x3,x4,x5,x6,x7,x8]|geometric=1);
                   1087:
1.6     ! takayama 1088: // 答.
1.1       takayama 1089: [14449864949304/9556267369631,10262588586540/9556267369631,13512615942680/9556267369631,
                   1090:  81112808747006/9556267369631,21816297744346/9556267369631,30858636683482/9556267369631,
                   1091:                               25258717886900/9556267369631,51191421070148/9556267369631]
                   1092:
                   1093:
                   1094: /*
1.6     ! takayama 1095:  上のデータで 0 を 1 に変更.
1.1       takayama 1096:   2 1 1
                   1097:   8 3 3
                   1098:   1 2 6
                   1099:
                   1100:   row sum: 4,14,9
                   1101:   column sum: 11,6,10
                   1102: */
1.6     ! takayama 1103: // 3 x 3 分割表.
1.1       takayama 1104: A=[[0,0,0,1,1,1,0,0,0],
                   1105:    [0,0,0,0,0,0,1,1,1],
                   1106:    [1,0,0,1,0,0,1,0,0],
                   1107:    [0,1,0,0,1,0,0,1,0],
                   1108:    [0,0,1,0,0,1,0,0,1]];
                   1109: B=[14,9,11,6,10];
                   1110: hgm_ahg_expected_values_contiguity(A,B,[1,1/2,1/3,1,1/5,1/7,1,1,1],[x1,x2,x3,x4,x5,x6,x7,x8]|geometric=1);
                   1111:
1.6     ! takayama 1112: // 期待値, 答.
1.1       takayama 1113: [207017568232262040/147000422096729819,163140751505489940/147000422096729819,217843368649167296/147000422096729819,
                   1114:  1185482401011137878/147000422096729819,358095302885438604/147000422096729819,514428205457640984/147000422096729819,
                   1115:  224504673820628091/147000422096729819,360766478189450370/147000422096729819]
                   1116:
1.6     ! takayama 1117: // Z やその微分の計算は hgm_ahg_contiguity 関数がおこなうが, これの簡易インターフェースは
        !          1118: // まだ書いてない.
1.1       takayama 1119:
                   1120:
1.6     ! takayama 1121: 4. x_ij は [GM2016] の1章で,
        !          1122:  たとえば 3x3 の時 [[1,1,1],[x_11,x_12,1],[x_21,x_22,1]]
        !          1123: となっているが, [GM2016] の Prop 7.1 の対応では,
        !          1124:    p = [[1,x_11,x_12],[1,x_21,x_22],[1,1,1]] となっているので注意.

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>