Part 4: Higher Mathematics / より高度な数学

Matrices and linear algebra / 行列と線形代数

Matrices are represented by arrays of rank 2, whereas vectors (in the more general sense of linear algebra, not geometrical vectors) are represented by arrays of rank 1. The module Numeric contains a few operations that are particularly useful for dealing with matrices and vectors:

行列は、階数２の配列によって表現されます。一方、ベクトル（幾何ベクトルではなく、線形代数におけるより一般的な意味でのベクトル）は、階数１のベクトルで表現されます。 Numeric モジュールには、行列やベクトルを扱うのに実用的な演算がいくつかあります。

dot(a, b) returns the dot product of a and b. This can be used with two vectors (scalar product), a matrix and a vector (matrix-vector product), or two matrices (matrix product). The last dimension of a must have the same length as the first dimension of b.
transpose(a) returns the transpose of the matrix a. This function can actually be applied to arrays of any rank and in general reverses the order of the dimensions.
identity(n) returns an nxn unit matrix.
ravel(a) returns a vector containing all elements of the array a. The order of the elements is defined by the first index varying most slowly, i.e. the elements of a 2x2 array are returned in the order (0,0), (0,1), (1,0), (1,1).
reshape(a, shape) returns an array of the specified shape (a tuple) containing the elements of the array a (in the same order as in ravel(a)). This is often useful to construct special matrices. The shape must lead to the same total number of elements. This restriction is lifted in the function resize(a, shape), which allows any shape. If more elements are needed than are contained in a, they are repeated from the start.
dot(a, b) は、a と b のドット積を返します。これは、二つのベクトルのスカラー積、行列とベクトルの積、二つの行列の積、に使えます。 a の最後の次元と b の最初の次元の長さは等しくなくてはなりません。
transpose(a) は、行列 a の転置行列を返します。実際には、この関数は任意の階数の配列に適用でき、一般的には次元の順序を反転させます。
identity(n) は、n 行 n 列の単位行列を返します。
ravel(a) は、配列 a の全ての要素からなるベクトルを返します。要素の順序は、最初の添字が最も緩やかに変化するように決められます。すなわち、２行２列の配列要素は、(0,0), (0,1), (1,0), (1,1) の順序で返されます。
reshape(a, shape) は、配列 a の要素からなり、（タプルで）指定された形状を持つような配列を返します。（配列要素の順序は ravel(a) と同じです。）これは、特殊な行列を作る際にしばしば便利です。この制限は resize(a, shape) 関数では任意の形状へと緩和されています。a が持つよりも多くの要素が必要な場合には、最初の要素に戻って繰り返されます。

The module LinearAlgebra contains various common operations on matrices. Matrices are represented by two-dimensional arrays. Some operations require square matrices, and check their input accordingly. All operations are available for real and for complex matrices. The functions are based on LAPACK-routines, which are well tested and efficient.

LinearAlgebra モジュールには、様々な慣用の行列演算があります。行列は二次元の配列で表されます。正方行列でなければならないために入力のチェックがなされるような演算もあります。全ての演算は、実行列にも複素行列にも使えます。これらの関数は、十分にテストされ効率の良い LAPACK ルーチンを基礎に作られています。

The function solve_linear_equations(a, b) solves a system of linear equations with a square non-singular matrix a and a right-hand-side vector b. Several right-hand-side vectors can be treated simultaneously by making b a two-dimensional array (i.e. a sequence of vectors). The function inverse(a) calculates the inverse of the square non-singular matrix a by calling solve_linear_equations(a, b) with a suitable b.

solve_linear_equations(a, b) 関数は、正方正則行列 a と右側ベクトル b からなる線形方程式系の解きます。 b を二次元配列（すなわちベクトルの列）とすることによって、複数の右ベクトルを同時に扱うことができます。 inverse(a) 関数は、適切な b を用いて solve_linear_equations(a, b) を呼び出すことにより、正方正則行列 a の逆行列を計算します。

The function eigenvalues(a) returns the eigenvalues of the square matrix a. The function eigenvectors(a) returns both the eigenvalues and the eigenvectors, the latter as a two-dimensional array (i.e. a sequence of vectors).

eigenvalues(a) 関数は、正方行列 a の固有値を返します。 eigenvectors(a) 関数は、固有値と固有ベクトルを返し、後者は二次元配列（すなわち、ベクトルの列）です。

The function determinant(a) returns the determinant of the square matrix a.

determinant(a) 関数は、正方行列 a の行列式を返します。

The function singular_value_decomposition(a) returns three arrays V, S, and W^T whose matrix product is the original matrix a. V and W^T are unitary matrices (rank-2 arrays), whereas S is the vector (rank-1 array) of diagonal elements of the singular-value matrix. This function is mainly used to check whether (and in what way) a matrix is ill-conditioned.

singular_value_decomposition(a) 関数は、 V, S, W^T の三つの配列を返します。これらの行列積は元の行列 a です。 V と W^T はユニタリ行列（階数２の配列）であり、一方 S は特異値行列の対角要素からなるベクトル（階数１の配列）です。この関数は、行列が悪条件にないかどうか（あるとしたらどのようにか）をチェックするのに主に使われます。

The function generalized_inverse(a) returns the generalized inverse (also known as pseudo-inverse or Moore-Penrose-inverse) of the matrix a. It has numerous applications related to linear equations and least-squares problems.

generalized_inverse(a) 関数は、一般化逆行列（偽逆行列またはムーア・ペンローズ逆行列としても知られています）を返します。これには、線形方程式や最小二乗法の問題に関連した応用が沢山あります。

The function linear_least_squares(a, b) returns the least-squares solution of an overdetermined system of linear equations. An optional third argument indicates the cutoff for the range of singular values (defaults to 10^-10). There are four return values: the least-squares solution itself, the sum of the squared residuals (i.e. the quantity minimized by the solution), the rank of the matrix a, and the singular values of a in descending order.

linear_least_squares(a, b) 関数は、データ過剰な線形方程式系の最小二乗解を返します。第三の引数はオプションで、特異値の範囲を決めるカットオフを与えます（デフォルトは 10^-10）。返り値は四つあります：降順で、最小二乗解そのもの、残渣の二乗和（すなわち、解が最小にする量）、行列 a の階数、そして a の特異値です。

Example / 例：

import Numeric; N = Numeric
import LinearAlgebra; LA = LinearAlgebra

# Create a matrix
# 行列の生成
a = N.reshape(N.arrayrange(25.), (5, 5)) + N.identity(5)

# Calculate the inverse
# 逆行列の計算
inv_a = LA.inverse(a)

# Verify the inverse by printing the largest absolute element
# of a * a^{-1} - identity(5)
# a * a^{-1} - identity(5) の最大絶対値の要素を印字し、逆行列を検証
print "Inversion error:", \
      N.maximum.reduce(N.fabs(N.ravel(N.dot(a, inv_a)-N.identity(5))))

# Print the eigenvalues of the matrix
# 行列の固有値を印字
print "Eigenvalues: ", LA.eigenvalues(a)

Note the way the abbreviations for the module names are created. Modules are just special data objects, so it is possible to assign them to variables.

モジュール名の略記名を作る方法に注意して下さい。モジュールは特殊なデータオブジェクトに過ぎないので、変数に割り付けることも可能なのです。

Fourier Transforms / フーリエ変換

The module FFT contains functions that perform various variants of the Fast Fourier Transform in one and two dimensions, using the package FFTPACK. All functions can be applied to arrays of any rank; the dimension(s) being transformed can be specified and default to the last one or two.

FFT モジュールには、一および二次元の種々の高速フーリエ変換を行う関数があります。FFTPACK というパッケージを使っています。全ての関数は任意の階数の配列に適用できます；変換すべき次元のデフォルトは最後の一つあるいは二つですが、それ（ら）を指定することも可能です。

The function fft(a) performs a standard complex one-dimensional FFT, and the function inverse_fft(a) performs the inverse transform. Both functions take an optional second argument specifying the length of the data sequence (by default the length of the array); if it is larger than the array, an appropriate number of zeros will be added. An optional third argument specifies the dimension along which the transform will be performed.

関数 fft(a) は、標準的な複素１次元 FFT を、inverse_fft(a) は逆変換を行います。両者とも、オプションでデータ列の長さ（デフォルトは配列の長さ）を指定する第二引数を受け取ります；それが配列長よりも大きい場合、適当数のゼロを追加します。オプションの第三引数は、変換を実行する際の変数となるべき次元を指定します。

The function fft2d(a) performs a two-dimensional FFT. An optional second argument specifies the shape (a tuple of length two) of the subarray to be transformed; by default the whole array is transformed. An optional third argument (also a tuple of length two) specifies the two dimensions for the transformation.

関数 fft2d(a) は、二次元 FFT を行います。オプションの第二引数により、変換すべき部分配列の形状 (shape)（長さ２のタプル）を指定します；デフォルトでは、配列全体が変換されます。オプションの第三引数（これも長さ２のタプル）は、変換の二変数の次元を指定します。

Parameter fitting / パラメータ最適化

The module LeastSquares provides the function leastSquaresFit(model, parameter, data) that performs general non-linear least-squares fit using the Levenberg-Marquardt algorithm.

LeastSquares モジュールには、関数 leastSquaresFit(model, parameter, data) が提供されています。これは、Levenberg-Marqardt アルゴリズムを用いた一般的な非線形最小二乗法最適化を行います。

The argument model specifies the function to be fitted. It will be called with two arguments: the first is a tuple containing the fit parameters, and the second is the first element of a data point (see below). The return value must be a number.

引数 model は、当てはめるべき関数を指定します。これは二つの引数を渡して呼び出します：一番目は最適化すべきパラメータからなるタプル、二番目はデータ点の第一要素（下記参照）です。返り値は一つの数です。

The argument parameter is a tuple of initial values for the fit parameters. The speed of convergence of the method depends on the quality of the initial values, so they should be chosen with care.

引数 parameter は、最適化するパラメータの初期値からなるタプルです。この方法の収束速度は、初期値の性質に依存するので、注意して選ばなくてはなりません。

The argument data is a sequence of data points to which the model is to be fitted. Each data point is a sequence of length two or three. Its first element specifies the independent variables of the model. It is passed to the model function as its first parameter, but not used in any other way. The second element of each data point tuple is the number that the return value of the model function is supposed to match as well as possible. The third element (which defaults to 1.) is the statistical variance of the data point, i.e. the inverse of its statistical weight in the fitting procedure.

引数 data は、モデル関数を適合すべきデータ点の列です。各々のデータ点は長さ２または３の列からなります。その第一要素は、モデルの独立変数を指定するもので、モデル関数に第一パラメータとして渡され、それ以外の用途には使われません。各データ点のタプルの第二要素は、モデル関数の返値が出来るだけ一致すべき対象となる数値です。第三要素（デフォルト値は１）はデータ点の統計的分散、すなわち、最適化の手続きにおけるその点の統計的重みの逆数です。

Example / 例：

from LeastSquares import leastSquaresFit
import Numeric

# The mathematical model.
# 数学的モデル
def exponential(parameters, x):
    a = parameters[0]
    b = parameters[1]
    return a*Numeric.exp(-b/x)

# The data to be fitted to.
# 適合化されるべきデータ
data = [(100, 4.999e-8),
	(200, 5.307e+2),
	(300, 1.289e+6),
	(400, 6.559e+7)]

fit = leastSquaresFit(exponential, (1e13, 4700), data)

print "Fitted parameters:", fit[0]
print "Fit error:", fit[1]

Functions on a grid / 格子点で与えられた関数

In numerical calculations many functions are defined by a values on a grid of evaluation points. Evaluation of such functions for off-grid points requires interpolation. The module Scientific.Functions.Interpolation defines such functions and provides some common operations on them. The grids must be orthogonal, but do not have to be equidistant. A grid is defined by a set of axes, one per dimension, specifying the coordinates of the grid points.

数値計算では、関数の多くは格子点上で計算された一組の値で定義されます。格子点以外でこのような関数の値を求めるには、内挿が必要です。 Scientific.Functions.Interpolation というモジュールは、このような関数を定義し、よく使われる演算を提供します。格子は直交でなくてはなりませんが、等間隔である必要はありません。格子は各点の座標を指定する一組の軸で定義されます。一次元あたり軸が一つです。

A function is defined by InterpolatingFunction(grid, value), where grid is a sequence of axes whose length determines the dimensionality of the grid. The argument values must be an array with a rank equal to or higher than the dimensionality of the grid. If it is higher, then the function has multiple values (e.g. a vector field). An optional third argument specifies a value for the function outside the grid. If none is given, evaluation of the function outside the grid results in an error.

関数は、InterpolatingFunction(grid, value) で定義します。ここで、grid は軸からなる列で、その長さによって格子の次元が決まります。引数 values は、格子の次元と同じまたはより高い階数を持つ配列です。階数がより高い場合には、関数は多価です（例えば、ベクトル場）。第三引数はオプションで、格子の外側の関数値を指定します。なにも与えられない場合、格子の外側で関数値を求めるとエラーになります。

An interpolating function can be called like any other function. The number of arguments must be equal to the dimensionality of the grid, and the arguments must be integers or real numbers.

内挿関数は、他の任意の関数と同様に呼び出すことができます。引数の個数は、格子の次元と等しくなくてはなりません。また、引数の型は整数または実数です。

Interpolating functions offer some useful operations in addition to a simple evaluation:

内挿関数には、単純な求値だけでなく、幾つかの便利な演算操作が提供されています：

function.derivative(variable) returns another interpolating function that is the (numerical) derivative with respect to the variable specified (default to zero).
function.integral(variable) returns another interpolating function that is the (numerical) indefinite integral with respect to the variable specified (default to zero).
function.definiteIntegral(variable) returns another interpolating function that is the (numerical) definite integral with respect to the variable specified (default to zero). This operation reduces the dimensionality by one; for a one-dimensional function, a simple number is returned.
function.selectInterval(first, last, variable) returns an interpolating function defined on a smaller grid, reduced to the interval from first to last along dimension variable. This is most useful for the calculation of definite integrals.
function.derivative(variable) は、variable 番目の変数（デフォルトはゼロ）に関する（数値的）微分を表す別の内挿関数を返します。
function.integral(variable) は、variable 番目の変数（デフォルトはゼロ）に関する（数値的）不定積分を表す別の内挿関数を返します。
function.definiteIntegral(variable) は、variable 番目の変数（デフォルトはゼロ）に関する（数値的）定積分を表す別の内挿関数を返します。この操作は、次元を一つ減らします；一次元関数の場合は、一個の数値が返されます。
function.selectInterval(first, last, variable) は、variable 番目の次元に沿って first から last までの間隔に縮小された、より小さな格子で定義された内挿関数を返します。これは、定積分の計算において最も有用です。

Example / 例：

from Scientific.Functions.Interpolation import InterpolatingFunction
import Numeric; N = Numeric

# Define an axis with 10 points
# 軸 (axis) を１０個の点で定義
axis = N.arrayrange(0, 1.1, 0.1)

# Calculate the values (here a simple square root)
# 値の計算（ここでは単純な平方根）
values = N.sqrt(axis)

# Create the interpolating function
# 内挿関数の生成
f = InterpolatingFunction((axis,), values)

# Print a table of values for comparison
# 比較のために値の表を印字
for x in N.arrayrange(0, 1., 0.13):
    print "%5.2f: %10.5f %10.5f" % (x, N.sqrt(x), f(x))

Visualization / 可視化 (VRML)

Simple plots are not always sufficient for visualizing complex data, especially in three dimensions. The module Scientific.Visualization.VRML provides a straightforward way to construct three-dimensional scenes composed of standard objects (cubes, spheres, etc.). These scenes can be written to files in the VRML (Virtual Reality Modelling Language) format, to be viewed with a VRML viewer. VRML viewers are available for most systems. For more information on VRML and VRML viewers, see the VRML Repository at the San Diego Supercomputer Center.

複雑なデータ（特に三次元）を描出するには、単純な図示では必ずしも十分ではありません。 Scientific.Visualization.VRML モジュールを用いれば、（立方体や球などの）標準的なオブジェクトからなる三次元的情景（シーン）を直ちに構築することができます。これらのシーンは、VRML (Virtual Reality Modelling Language) フォーマットでファイルに書き出し、VRML 表示ソフトで見ることができます。 VRML 表示ソフトは大抵のシステムで使用できます。 VRML と VRML 表示ソフトについての詳細は、サンディエゴ・スーパーコンピュータセンターの VRML Repository を参照してください。

If you want to display VRML scenes directly from your Python programs, instead of writing a file and running the VRML viewer manually, you must set the Unix environment variable VRMLVIEWER to the name of the VRML viewer you have installed.

一々ファイルを書き出して VRML 表示ソフトを起動する代りに、Python プログラムから直接 VRML シーンを表示したければ、インストールした VRML 表示ソフトの名前を Unix の環境変数 VRMLVIEWER に設定する必要があります。

Note: this course covers only the basic use of the modules VRML. For more details see the source code.

注：このコースでは、VRML モジュールの基本的使用法のみ扱います。詳細については、ソースコードを見てください。

The typical sequence of operations is the following:

Create a VRML scene with Scene(objects), where objects is a list of VRML objects. You can pass an empty list to create an empty scene.
Create VRML objects (see below) and add them to the scene with scene.addObject(object).
Write the scene to a file with scene.writeToFile(filename), or view it directly using scene.view().

典型的な一連の操作は次のようになります：

VRML シーンを Scene(objects) によって生成。ここで、objects は、VRML オブジェクトのリスト。空のリストを渡して、空のシーンを生成することも可能。
VRML オブジェクトを生成（下参照）し、scene.addObject(object) でシーンに追加。
scene.writeToFile(filename) でシーンをファイルに書き出す。または、直接 scene.view() によって描画。

VRML objects are defined by geometrical data, e.g. a sphere is defined by a center (a vector) and a radius (a number). In addition, objects can have one or more properties that affect their appearance. The most important property is the material, which defines color, transparency, reflectivity, etc. In this course, only diffuse materials in various predefined colors are used, i.e. non-transparent non-reflecting materials. See a book on VRML for the full possibilities.

VRML オブジェクトは、幾何的データによって定義されます。例えば、球は中心を表す一つのベクトルと半径を表す一つの数で定義されます。さらに、オブジェクトはその外観に係わる一つ以上の性質 (properties) を持つことができます。最も重要な性質は、素材 (material) で、色、透明度、反射度などを定義します。このコースでは、予め決められた種々の色を持つ広汎的素材 (diffuse materials) 、すなわち、透明でもなく反射もない素材のみを使います。全ての可能な性質に関しては、VRML に関する書籍を参照して下さい。

Attributes are always optional and specified by keyword parameters. The material, for example, is specified by material=the_material. An object with no material specification has a white diffuse material by default.

属性は常にオプションで、キーワードパラメータによって指定します。例えば、素材は material=the_material のように指定します。指定の無い場合のデフォルトは、白の広汎的素材になります。

A material is created with DiffuseMaterial(color), where color is the color name. The available colors are black, white, grey, red, green, blue, yellow, magenta, and cyan. All of these can be prefixed with 'light ' or 'dark ' (note the space), e.g. 'light blue'.

素材は、DiffuseMaterial(color) によって生成されます。ここで、color は色の名前です。使用できる色は、black, white, grey, red, green, blue, yellow, magenta, cyan です。例えば 'light blue' のように、これら全色の先頭に 'light ' または 'dark '（空白に注意）を付けることができます。

The available VRML objects are:

使用できる VRML オブジェクトは以下の通りです：

Sphere(center, radius) creates a sphere.
Cube(center, edge) creates a cube centered around a given point. The edges are oriented along the x, y, and z axes.
Cylinder(point1, point2, radius) creates a cylinder whose axis runs from point1 to point2 with the radius specified.
Cone(point1, point2, radius) creates a cone whose base is a circle around point1 with radius radius and whose tip is at point2.
Line(point1, point2) creates a line from point1 to point2.
Arrow(point1, point2, radius) creates an arrow from point to point2 with a shaft radius of radius.
Sphere(center, radius) は、球を生成します。
Cube(center, edge) は、与えられた点を中心とする立方体を生成します。各辺は、x, y, z 座標に沿っています。
Cylinder(point1, point2, radius) は、軸が point1 から point2 までで、指定された半径を持つ円柱を生成します。
Cone(point1, point2, radius) は、point1 を中心とする半径 radius の円を底面とし、先端が point2 にある円錐を生成します。
Line(point1, point2) は、point1 から point2 までの線分を生成します。
Arrow(point1, point2, radius) は、軸半径が radius で、point1 から point2 を指す矢印を生成します。

The following example program reads a PDB file and constructs a scene consisting of a black sphere for each C-alpha atom and red lines connecting neighbouring C-alpha atoms.

次のプログラム例は、PDB ファイルを読み込み、各々が黒い球で表され、隣り合うものが赤い線で結ばれたα炭素からなるシーンを構築します。

from Scientific.IO.TextFile import TextFile
from Scientific.IO.FortranFormat import FortranFormat, FortranLine
from Scientific.Geometry import Vector
import string
from Scientific.Visualization import VRML

generic_format = FortranFormat('A6')
atom_format = FortranFormat('A6,I5,1X,A4,A1,A3,1X,A1,I4,A1,' +
                            '3X,3F8.3,2F6.2,7X,A4,2A2')

# Read the PDB file and make a list of all C-alpha positions
# PDB ファイルを読み込み、全てのα炭素原子の位置からなるリストを作る。
c_alpha_positions = []
for line in TextFile('protein.pdb'):
    record_type = FortranLine(line, generic_format)[0]
    if record_type == 'ATOM  ' or record_type == 'HETATM':
        data = FortranLine(line, atom_format)
        atom_name = string.strip(data[2])
        position = Vector(data[8:11])
	if atom_name == 'CA':
	    c_alpha_positions.append(position)

# Create the VRML scene and the materials and radii for spheres and lines
# VRML シーン、球の素材と半径、線の素材を生成。
scene = VRML.Scene([])
c_alpha_material = VRML.DiffuseMaterial('black')
c_alpha_radius = 0.2
line_material = VRML.DiffuseMaterial('red')

# Create the spheres and put them into the scene
# 球を生成し、シーンに配置。
for point in c_alpha_positions:
    scene.addObject(VRML.Sphere(point, c_alpha_radius,
				material=c_alpha_material))

# Create the lines and put them into the scene
# 線を生成し、シーンに配置。
for i in range(len(c_alpha_positions)-1):
    point1 = c_alpha_positions[i]
    point2 = c_alpha_positions[i+1]
    scene.addObject(VRML.Line(point1, point2, material=line_material))

# View the scene
# シーンを描画。
scene.view()

Other scientific modules / その他の科学計算モジュール

There are more modules for scientific applications which are not covered in this course. See the manual for documentation.

このコースで網羅されていない科学計算モジュールは、まだまだあります。マニュアルを参照して下さい。

Exercises / 練習問題

Write a function that applies a mathematical function (exp, sqrt, etc.) to a square matrix by diagonalization. The steps are: 1) find the eigenvalues and eigenvectors of the matrix 2) apply the mathematical function to the eigenvalues 3) undo the diagonalization by multiplying with the eigenvector matrix from both sides.
Test this program by calculating the exponential of a matrix alternatively from the Taylor expansion.
Write a program that reads two PDB files containing two conformations of the same protein and produces a VRML scene containing an arrow from conformation 1 to conformation 2 for each C-alpha atom.
Hint: put the code for extracting the C-alpha positions from a file into a function and call this function twice.
対角化法により、正方行列に算術関数 (exp, sqrt など) を適用する関数を書きなさい。手順は次の通りです： 1) 行列の固有値と固有ベクトルを見つける。 2) 固有値に算術関数を適用する。 3) 固有ベクトル行列を両側から掛け、対角化の逆操作を行う。
テーラー展開による別法で行列の指数関数を計算し、このプログラムをテストしなさい。
同一のタンパク質の異なる二配置を表す二つの PDB ファイルを読み込み、各α炭素原子について配置１から配置２へ向う矢印を描いた VRML シーンを制作するプログラムを書きなさい。
ヒント：ファイルからα炭素の位置を取り出すプログラムを関数化し、これを二回呼ぶようにしましょう。

Table of Contents / 目次