/*************************************************************/
/*              コレスポンデンス分析 : CORRESP               */
/*                                                           */
/* DATE: 1997/07/27                                          */
/* CODE: Tokuhisa SUZUKI                                     */
/*************************************************************/


proc iml ;

/*-------------------------------------------------------------
                   データ行列 z の指定
-------------------------------------------------------------*/

    z = { 1 1 2 ,
          0 3 1 ,
          2 2 4 ,
          1 0 3
        }
    ;

/*-----------------------------------------------------------*/

    p  = nrow( z )         ;  /* 行の要素の個数 */
    q  = ncol( z )         ;  /* 列の要素の個数 */
    r  = diag( z[,+] )     ;  /* 行和の対角行列 */
    c  = diag( z[+,] )     ;  /* 列和の対角行列 */
    rs = inv( sqrt( r ) )  ;
    cs = inv( sqrt( c ) )  ;
    x  = rs * z * cs       ;  /* Ｚを規準化した行列Ｘ */
    call svd( u, l, v, x ) ;  /* Ｘの特異値分解       */
    a = rs * u * diag( l ) ;  /* サンプル・スコア行列 */
    b = cs * v * diag( l ) ;  /* カテゴリ・スコア行列 */

    reset noname ;
    print / l[ format = 10.6  colname = "特異値"           ] ;
    print / a[ format = 10.3  colname = "サンプル・スコア" ] ,
            b[ format = 10.3  colname = "カテゴリ・スコア" ] ;

    out = a[ , 2:3 ] // b[ , 2:3 ] ;   /* トリビアルな解を除去 */
    _type_ = repeat( {"OBS"}, p, 1 ) // repeat( {"VAR"}, q, 1 ) ;
    create coord from out[ rowname = _type_ ] ;
    append       from out[ rowname = _type_ ] ;
    close  coord ;
quit ;

data coord ;
   keep col1 col2 _type_ _name_  ;
   length _type_ $3 _name_ $16;
   rename col1 = dim1 col2 = dim2 ;
   set coord ;
   select( _type_ ) ;
       when( "OBS" ) do ;
           n + 1 ;
           _name_ = left( put( n, 3. ) ) ;     /* 自動サンプル番号 */
       end ;
       when( "VAR" ) do ;
          m + 1 ;
          _name_ = "V" || left( put( m, 3. ) ) ; /* 自動変数ラベル */
       end ;
   end ;
run ;

data anno ;
   set coord ;
   length function $8 text $16;
   xsys  = '2'        ;
   ysys  = '2'        ;
   text  = _name_     ;
   style = 'HWDMX003' ; /* ＭＳ明朝 = HWDMX003  */
   x     = dim1       ;
   y     = dim2       ;

   if _type_ = 'OBS' then
       do;
           size     = 1.0     ; /* サンプル番号の印刷サイズ */
           color    = 'BLACK' ;
           position = '5'     ;
           function = 'LABEL' ;
           output ;
       end;

   if _type_ = 'VAR' then
       do;
           color = 'BLACK ';
           size  = 1.0     ;  /* 変数名の印刷サイズ */
           if dim1 >= 0 then position = '6' ; /* left justify  */
           else              position = '4' ; /* right justify */
           function = 'LABEL';
           output ;
       end;
run ;

goptions hsize   = 17 cm  vsize   = 20 cm
         horigin =  2 cm  vorigin =  3 cm
         device = winprtm
;
proc gplot data = coord ;
    plot dim2 * dim1 / anno = anno frame
         href   = 0     vref   = 0
         lvref  = 3     lhref  = 3
         vaxis  = axis2 haxis  = axis1
         vminor = 1     hminor = 1
    ;
    axis1  length = 13 cm
           offset = ( 2       )
           value  = ( h = 0.5 )
           label  = ( h = 0.5 )
           order  = ( -0.15 to 0.25 by 0.05 ) /* データに依存 */
    ;
    axis2  length = 13 cm
           offset = ( 2 )
           value  = ( h = 0.5 )
           label  = ( h = 0.5 )
           order  = ( -0.15 to 0.25 by 0.05 ) /* データに依存 */
    ;
   symbol v = none ;
   footnote1 f = hwdmx003 h = 0.7 '図ｘ　コレスポンデンス分析' ;
run  ;
quit ;
footnote1 ;