ここでは、JSLを使って作成できる特殊な行列（単位行列や対角行列など）について説明します。
Identity()関数を使うと、任意の次元で単位行列を作成できます。単位行列とは、対角線上の要素が1で、それ以外の要素が0の正方行列です。引数には次元数のみを指定します。
Identity( 3 );
[1 0 0,
0 1 0,
0 0 1]
関数J()は、第1引数を行数、第2引数を列数、第3引数を全要素の値とする行列を作成します。
J( 3, 4, 5 );
[5 5 5 5,
5 5 5 5,
5 5 5 5]
J( 3, 4, Random Normal() ); // 結果は異なる
[0.407709113182904 1.67359154091978 1.00412665221308 0.240885679837327,
-0.557848036549455 -0.620833861982722 0.877166783247633 1.50413740148892,
-2.09920574748608 -0.154797501010655 0.0463943433032137 0.064041826393316]
Diag()関数は、（同数の行と列を持つ）正方行列やベクトルから対角行列を作成します。対角行列とは、非対角要素がすべて0の正方行列です。
D = [1 -1 1];
Diag( D );
[1 0 0,
0 -1 0,
0 0 1]
Diag([1, 2, 3, 4]);
[1 0 0 0,
0 2 0 0,
0 0 3 0,
0 0 0 4]
A = [1 2,3 4];
f = [5];
D = Diag( A, f );
[1 2 0
3 4 0
0 0 5]
三つ目の例では、一見すると、すべての非対角要素がゼロというわけではありません。行列表記を使用すると、次のようになります。
[A 0,
0` f]
ここで、Aとfは例の行列で、0はゼロの列ベクトルです。
VecDiag()関数は、行列の対角要素から成るベクトルを作成します。
v = Vec Diag(
[1 0 0 1, 5 3 7 1, 9 9 8 8, 1 5 4 3]);
[1, 3, 8, 3]
Vec Quadratic()関数は、回帰分析における予測の標準誤差や、外れ値に関するMahalanobis距離やT2乗統計量などを計算します。Vec Quadratic(Sym, X)は、Vec Diag(X*Sym*X`)を計算するのと同じです。第1引数は対称行列でなければなりません。通常は、共分散行列の逆行列を指定します。第2引数は、列数が第1引数の行列と等しい矩形行列です。
Trace()関数は、正方行列の対角要素の和を戻します。
D = [0 1 2, 2 1 0, 1 2 0];
Trace( D ); // 1を戻す
Index()関数は、最初の引数から最後の引数までの整数値を要素とする行ベクトルを作成します。2重コロン::は等価の二項演算子です。
6::10;
[6 7 8 9 10]
Index( 1, 5 );
[1 2 3 4 5]
オプションの第3引数を指定して、増分をデフォルト値の+1から変更することもできます。
Index( 0.1, 0.4, 0.1 );
[0.1, 0.2, 0.3, 0.4]
増分は負の数でもかまいません。
Index( 6, 0, -2 );
[6, 4, 2, 0]
デフォルトの増分は1、または、最初の引数が2番目の引数よりも大きい場合は-1です。
Shape()関数は、既存の行列を指定の次元に作り直します。たとえば次のように指定すると、3x4の行列であるaが12x1の行列に変わります。
a = [1 1 1, 2 2 2, 3 3 3, 4 4 4];
Shape( a, 12, 1 );
[1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4]
行列における各列の要約統計量を、行ベクトルで戻す関数がいくつか用意されています。
mymatrix = [11 22, 33 44, 55 66];
V Max( mymatrix ); // 各列の最大値を戻す
[55 66]
V Min( mymatrix ); // 各列の最小値を戻す
[11 22]
V Mean( mymatrix ); // 各列の平均値を戻す
[33 44]
V Median( mymatrix ); // 各列の中央値を戻す
[33 44]
V Quantile( mymatrix, 0.75 ); // 各列の75%点を戻す
[55 66]
V Sum( mymatrix ); // 各列の合計を戻す
[99 132]
V Std( mymatrix ); // 各列の標準偏差を戻す
[22 22]
Design()関数は、ベクトルまたはリストの計画行列を作成します。引数の一意の値ごとに1つの列を作成します。計画行列の要素は0または1です。たとえば、下のxは、1、2、3の値から成っているため、計画行列には1の列、2の列、そして3の列があることになります。行列の各行には、行の値に対応する列に1が入ります。次の例では、第1行（1）に対しては、１の列（第1列）に1が入り、他の列には0が入ります。第2行（2）に対しては、2の列に1が入り、他の列には0が入ります。以下、同じような操作で行列が作成されます。
x = [1, 2, 3, 2, 1];
Design( x );
[1 0 0,
0 1 0,
0 0 1,
0 1 0,
1 0 0]
この変形としてDesign Nom()関数またはDesign F()関数があり、これは最後の列を削除して、それを他の列の値から引きます。したがって、Design Nom()またはDesign F()でできる行列の要素は、0、1、または–1となります。また、Design Nom()またはDesign F()でできた行列の列の数は、元のベクトルが持つ要素の種類の数より1つ少なくなります。この関数は、効果のフルランクの計画行列を作成します。
x = [1, 2, 3, 2, 1];
Design Nom( x );
[1 0, 0 1, -1 -1, 0 1, 1 0]
Design Nom()については、分散分析（ANOVA）の例でさらに説明しています。
順序尺度の因子を簡単にコーディングするには、Design Ord()関数を使います。この関数は、ベクトルまたはリストの一意の値の数より1つ少ない列を使ってフルランクのコーディングを作成します。低水準の行をすべてゼロにし、それに続く水準は計画行列の行にひとつずつ1を加えていきます。
x = [1, 2, 3, 4, 5, 6];
Design Ord( x );
[0 0 0 0 0,
1 0 0 0 0,
1 1 0 0 0,
1 1 1 0 0,
1 1 1 1 0,
1 1 1 1 1]
Design()、Design Nom()、およびDesign Ord()では、全水準の値とその順序を指定する第2引数を取ることができます。この機能により、たった1行に対する計画行列も作成することができます。
• Design( values, levels )は計画行列を作成します。
• Design Nom( values, levels )はフルランクの計画行列を作成します。
メモ:
• values引数は、単一の要素でも、要素の行列またはリストでもかまいません。
• levels引数は、計画行列作成の対象となる全水準を値として持つリストまたは行列です。
• 作成される行列は、values引数として指定した値の要素数と同じ行数になります。
• 作成される行列は、levels引数として指定した項目数と同じ列数になります。Design Nom()とDesign Ord()の場合は、levels引数として指定した項目数より1つ少ない列数になります。
• 値が見つからない場合は、行全体がゼロになります。
Design( 20, [10 20 30] );
[0 1 0]
Design( 30, [10 20 30] );
[0 0 1]
Design Nom( 20, [10 20 30] );
[0 1]
Design Nom( 30, [10 20 30] );
[-1 -1]
Design Ord( 20, [10 20 30] );
[1 0]
Design( [20, 10, 30, 20], [10 20 30] );
[0 1 0,
1 0 0,
0 0 1,
0 1 0]
Design({"b", "a", "c", "b"}, {"a", "b", "c"});
[0 1 0,
1 0 0,
0 0 1,
0 1 0]
Direct Product()関数は、2つの正方行列の直積（Kronecker積）を求めます。行列と行列の直積は行列となり、その要素は、AとBの要素の積となります。
G = [1 2, 3 5];
H = [2 3, 5 7];
Direct Product( G, H );
[2 3 4 6,
5 7 10 14,
6 9 10 15,
15 21 25 35]
H Direct Product()関数は、行数の等しい2つの行列の行ごとの直積を求めます。
H Direct Product( G, H );
[2 3 4 6, 15 21 25 35]
HDirect Product()は、計画行列における交互作用の列を作る際に便利です。
XA = Design Nom( A );
XB = Design Nom( B );
XAB = HDirect Product( XA, XB );
X = J( NRow( A ), 1 ) || XA || XB || XAB;