つれづれ思うこと Perl 編 その８

では次に、最近少しは慣れてきた Mathematica について書いてみたいと思う。おっとまずい、どんどん Perl の話題から遠ざかっている。両者の間に何か関連はないだろうか？

このソフトは高価なのですが、無料で楽しめる方法が少なくとも２つあります。一つは、Raspberry pi です。これに Mathematica を無料でインストールできるのです。その代わり速度は遅く、10 倍ぐらいの処理時間がかかってしまうかもしれません。Mathematica は購入しても一つの PC にだけしかインストールできない仕組みになっています。そのため、筆者も家では Raspberry pi にのんびりと計算をさせ、翌日に職場の PC で本格的に計算させていました。Raspberry pi の OS には 32-bit と 64-bit があり、Mathematica は前者にだけインストール可能でした。しかし、今 HP を見てみると、2023/May からは 64-bit でも使えるようなことが書かれています。64-bit 版ですと、メモリーを 8 GB などと積めるので、Mathematica の動きも速くなるかもしれません。まだ試したことがないので、時間の空いた時に試してみたいと思います。

もうひとつ無料で楽しめる方法は、Jupyter lab を使う方法です。実は Mathematica という名前が指しているのは GUI 部分（という表現は変ですが）のことで、実際に計算しているコア部分は Wolfram-engine と呼ばれています。その engine 部分が無料で PC にインストールできるようになりました。すると、Windows のコマンドプロンプトなどから、計算を打ち込むことができます。しかし、いくらなんでも、これは使い勝手が悪い。そこで engine を Jupyter lab につなげると、ちょうど Python を Jupyter lab で操作するのと同じように Wolfram-engine を使うことができるようになります。もちろん Mathematica の方がインターフェースとしてはよいのですが、この Jupyter lab は Mathematica とよく似ているので（というより、Mathematica を見本として作られた？）まずまずのことができます。なお、計算能力は全く同じです。Mathematica ではできて、Wolfram-engine & Jupyter lab ではできないということはありません。今、家では Raspberry pi をやめて、この第二の方法を活用しています。Jupyter lab の上では、グラフをマウスで拡大したり回したりはできないのですが、これもコマンドベースで命令すればできますので、まあなんとかなります。PDF に出力させれば、論文用の高分解能の図も作れます。

この Mathematica は Lisp と呼ばれる言語で作られているらしく、いわゆる関数型言語に分類されるそうです。関数型が何かという専門的な部分は、私にはよく分かりません。使い始めてちょっと驚いたのは、例えば次のような IF 文です。

y = If[x < 0, -x, x]

これを普通？の言語で書くと、

if ($x < 0)
{
    $y = -$x;

}
else
{
    $y = $x;

}

となります（むりやり Perl にしてしまった）。。。

これだけだと特に変わった言語のようにも見えないのですが、もし、else 文の中に複数の行があると、次のようになります。「もし x が負の数ならば -x を y に代入する。逆に正の数ならば、x に 3 をかけて平方根をとり、これを y に代入する。」

y = If[x < 0, -x, z=x*3; Sqrt[z]]

よって、"z=x*3; Sqrt[z]" の部分がまとまって、まるで If 文の第三引数であるかのように振る舞うのです。この時の「;」が重要で、これで "z=x*3" と "Sqrt[z]" の二つが繋がって、あたかも一つの引数であるかのように振る舞います。よって、"Sqrt[z]" の後にまで「;」を付けてはいけません。もし付けてしまうと、もう一つ後ろに何らかの表現が続くことになってしまいます。

さらに、If 文が最終的に Sqrt[z] を返してくる、そして、それが y に代入されるという点も興味深いところです。このように If 文そのものが関数となっているわけです。Excel の =If もそうなっていたような気がする。

筆者は、どうも C, Perl などの癖が出てしまい、Mathematica でも以下のように書いてしまいます。各文末の「,」や「;」に気を付けないといけませんが、これですと一般的な手続き型言語に似た構文となります。ここでも Sqrt[z] の後に「;」を付けないように気をつけないといけません。

(* の箇所はコメントとして無視されます。 *)

If

[
    (* if *) x < 0,
        y = -x,
    (* else *)
        z = x*3;
        y = Sqrt[z]
]

このように If も For も全てが関数として作られており、この関数の引数の中に別の、あるいは同じ関数を入れ込むことができます。この「引数の中に関数を書く」という操作をどんどん繰り返していくと、入れ子（ロシアのマトリョーシカ人形）のようになります。そんな構文を何に使うのかというと、漸化式です。Mathematica においても For 文を使わずに、この漸化式（数学的帰納法）のコンセプトで書く方法があるそうですが、私にはとてもできません。筆者の知っている人で OCaml が大好きという方がおり、何故かと尋ねたところ、一つの理由はこの関数型にあるのだそうです。

つれづれ思うこと Perl 編 その７

まだ Python の NumPy の存在に気づかなかった 2000 年代前半の頃、Octave という名の数値計算用ソフトがあることを知り使い始めました。これは便利で、今まで Perl などを使って書いていた面倒な数値計算が、たいへん簡単に行えるようになりました。Octave では「for 文を使ったら負け」という表現がよく出てきます。つまり、行列の計算などをする際に、できるだけライブラリの関数を使うようにせよということです。そうすると高速で演算できるのですが、それを安易にうっかり for 文を使って行と列の要素をひとつずつ繰り返しながら計算しようものなら、たいへん遅くなってしまいますよという意味です。よって、たえず「この計算は行列の演算に置き換えることはできないだろうか？」と考えながらプログラミングする癖がつきます。この考え方は、特に Octave に限らず、NumPy でも Mathematica でも同じですね。

NMR では数値シミュレーションだけでなく、数式処理も時々必要になってきます。3+4=7 のような数値計算ではなく、cos^2(wt)+sin^2(wt) = 1 のような計算が数式処理です。ω や sin, cos などを残したまま計算するのです。これが出来るのは限られた幾つかのソフトだけで、有名なのは Maple と Mathematica です。幸い、前にいた職場にはキャンパスライセンスがあったので、両者とも個人的には研究費を払わずに使うことができました。それで、Maple の HP を見てみると、モミジが描かれており、赤毛のアンを連想させるような雰囲気の中で、メープルシロップのたっぷりかかったモミジ型のお菓子（お土産で有名ですが、これおいしい！）の連想も手伝って、何も迷わずにこれの DVD を取り寄せました。

ちなみに、こうして自分の過去を思い返してみると、書籍の表紙に描かれている動物が可愛いかどうか？あるいは、連想させる食べ物がおいしいかどうか？でその後に使うプログラミング言語を決めてきたようなところがあります。ちゃんと情報教育を受けてきた人は、もっと理論的な考えでもってツールを選ぶものですが、生物という全く異なる分野の教育しか受けたことのない人間は案外、雰囲気というもので重要なことを決めてしまっているのかもしれません。もちろん Raspberry-pi 大好きです。そして Strawberry-Perl を愛用しています。これもフルーツ系の食べ物です。Red-Hat は嫌いで、ペンギンの描かれた Linux 系が好きです。ソフトに限らず、何か製品を開発していく人は、そのような雰囲気や感情で物事を決めていく人間もいるということを知り、ちょっとネーミングや宣伝を工夫してみると大当たりしてしまうかもしれません。もし、Raspberry-pi が Blue-Hat というような名前で、ギャング（やくざ）が青い帽子をかぶったような絵がイメージキャラクターだったとしたら、小中学校で教育用として人気を博すことはまずなかったでしょう。

数年間ほど楽しんで Maple を使っていたのですが、ある時、NMR の HSQC をシミュレーションしようとして動かしたところ、超複雑な答が 100 行ぐらい吐き出されてくるのです。理想的には Ix という二文字で終わりのはずなのですが。。。これは数式処理でしか起こらず、しかも全く同じ計算をしているのに 10 回に 1 回ぐらいしか出現しない。どうして変な答がたまに出るのかが分からず、また半年ほどあれやこれやと悩んだのですが、最終的に分かったことは、exp(i H t) の計算がうまく行っていないということでした。どこかで書いたと思うのですが、この H は行列で (i H t) の固有関数と固有値がピタッと数式処理で出れば exp(i H t) の計算は超簡単になります。しかし、この行列の対角化ができないと判断されれば、テーラー展開に切り替わるのです。しかし、テーラー展開の答はあくまで近似値にすぎません。HSQC の数式シミュレーションでは途中で多数の項が出てきますが、それらは位相回しによって最終的にはほとんどの項が打ち消され簡単な答だけが残るのです。まさに三角関数の加法定理などが式の簡略化に使われます。しかし、近似値どうしでは、この最後のキャンセルが完全な０にならないわけですね。そして、結果として多数の大小変な項が残ってしまうのです。

これを Maple の開発部門に伝えたところ、先方でも再現されバグであると認めてくれました（認められて嬉しいのやら、逆に計算できなくなって悲しいのやら）。それでバグが修正されるのを待っていたのですが、その間に筆者が今の職場に異動してしまいました。ところが残念なことに、ここはキャンパスライセンスがないのです。もしかしたら、もうバグが修正されているかもしれないのですが、25 万円を払って購入し、それを確かめるのも勇気がいり、もしかしたら開発部門から「バグを見つけてくれたので感謝料として、あなたには Maple パックをモミジ型お菓子とともに無料進呈いたします！」みたいな知らせが来ないかな～などと期待していました。しかし、こちらから不躾にそんな要求をするわけにもいかず、泣く泣く Maple を断念しました。それで、なんだか変な狐が描かれた Mathematica を買うことに。今だに Mathematica にはあまり慣れず、Maple が恋しいのですが、キャンパスライセンスがない大学では仕方がありません。頑張って Mathematica に慣れることにしました。

Maple はシンタックスが C 言語とちょっと似ているので、C, Perl, Python などを知っている人にはとっつきやすいと思います。GUI もきれいで、ずっと使い続けたいソフトの一つです。いずれにしても、Mathematica と Maple という二大ソフトが競争してくれるお陰で、ともに開発が進み、価格も（高いとは思いますが）25 万円で済んでいるのでしょう。そして、幸いにもキャンパスライセンスのある大学の人はどんどんこれを利用して、自分を教育すべきです。筆者は高校の時から数学があまり得意ではなく、その所為もあって数学とはもっとも無縁の学科を選びました。しかし、不得意であっても腐らずに勉強していれば、上記のような優秀なソフトが数学の不得意部分を補ってくれます。例えば、上記で exp(i H t) を書きましたが、本当に賢い人は、これを紙と鉛筆だけで対角化して計算し、HSQC などをさっとシミュレートしてしまいます。私にはそれがとても出来そうにないので、ソフトの計算力を借りるわけです。

つれづれ思うこと Perl 編 その６

前回の記事で C 言語について追加で書いたのは、実は Tcl/Tk との関連を書きたかったためです。Tcl/Tk を知ったのは、Perl, Python を知った直後でした。NMR で有名なフーリエ変換ソフト NMRPipe にもこれが使われていました（今も動いています）。Tcl 言語は Perl や Python と同じような（コンパイル不要の）スクリプト言語で、一方の Tk は GUI（グラフィカル・ユーザー・インターフェース）を作るためのツールキットです。それまでは、ちょっとした GUI を作るのも面倒だったのですが、この Tk のおかげでだいぶん楽になりました。その Tk が Tcl と相性がよいので両者が一緒に使われるようになったのだと思います。

Tcl 言語はひとことで言うと、たいへん可愛いに尽きます。Perl や Python がそれぞれ立派なリャマと大蛇だとすると、Tcl はコアラのような感じ。例のシリーズ本では、お猿さんの図が表紙になっていますが、猿のイメージとは大違い。ぱっと見た感じではシェルに見えます。ですので、ちょっとした計算（x = 5+3）でも、set x [expr 5+3] のような面倒な書き方になります。しかし、逆に c-shell も知らない筆者にとっては、これが新鮮で可愛く映ったのでした。

しかし、これだけなら Tcl/Tk を使おうとは思わなったのですが、Tcl にはそれ以上に面白い点がありました。Tcl は C 言語のコードをコンパイルする時に、そこにライブラリーとして含めることができるのです。つまり、以下のようになります。複雑な計算の箇所は高速化が必要ですので C 言語で書いておきます。そして、そのコードの入力と出力のパラメータを Tcl 言語と紐づけします。そして、Tcl ライブラリーを含めてコンパイルします。できあがった実行可能バイナリーファイルですが、入出力を Tcl 言語にすることができるのです。これはたいへん便利。いろいろなパラメータで同じような計算をし、結果をグラフに表して比べたい時、そのパラメータの入出力部分は Tcl でいかようにも変えることができるのです。ですので、二度目、三度目のコンパイルは不要。もし、C 言語だけでそれをしようとすると、「パラメータ a を入力してください」などと画面で促し、それを scanf で取り込むという面倒なことになります。もし、数百種類のパラメータでそれをしないといけないとすると、それでもう終わりです。

25~30 年ほど経った先日、当時のプログラムを使わないといけなくなり、バイナリーを動かそうとしてみましたがダメでした。やはりコンパイルし直さないと。しかし、コンパイル時にかなりの量の warning と error が出てしまいました。仕様がかなり変わっているのですね（引数に int, float などを明記しないとダメなど）。半日ほどかけて修正し、やっとほんの少しの warning を残した状態で make できました。残余双極子相互作用値と PDB から蛋白質分子の磁場配向テンソルを求めるプログラムなのですが、無事に計算できました。C 言語のアメーバを使っているのですが、初期値を５種類つくっての計算は 0.5 秒ぐらいでしょうか？この程度のフィッティングであれば、どのような初期値でスタートしても、いつも同じ最適解に収束します。

もう最近は C 言語から離れてしまったので、Tcl も使わなくなってしまいました。それに（間違えているかもしれませんが）いちどバージョンを上げるのに有料になった時があり、それを機会にやめてしまいました。手元に 2010 年に発刊された書籍「Tcl and the Tk Toolkit, 2nd ed」がありますので、また始めてみてもよいかもしれません。

つれづれ思うこと Perl 編 その５

C 言語についてもう少し。

Perl で exp(i H t) などの計算を試してはみたものの、やはり 1995 年頃は遅くて使い物にはならなかった。そこで科学計算にはやはり C 言語を使っていた。その時に役立った本が「Mumerical Recipes in C」（の翻訳版）。ここにはいろいろな計算のためのコードが載っていて、それをそのまま使ってかなりの数値計算をすることができる。特に重宝したのが最適化（フィッティング）の処理である。データに非線形的にフィッティングするための方法では Levenberg–Marquardt 法が有名であるが、コードを打ち込むのが大変そうだったので、何かもっと簡単な方法はないかな？とこの本で探してみると「アメーバ法」という方法が紹介されていた。正式には Nelder–Mead 法と呼ぶらしい。コードも数百行ぐらいで閉じているので、頑張って間違わずに打ち込めば、すぐに動く。この方法のよい点はモデル式の微分が要らない、不連続な箇所があってもよいなどで、まあとにかく制限がほとんど無いので、どんな状況にも適用できる。欠点は時間がかかる、間違えた解にトラップされやすいなど。よって、予想できる解がすでにある時には、できるだけそこからスタートするとよいです。

こんな（失礼ながら）子供遊びのような方法で数学的な最適化ができるのか？と思われ勝ちだが、試してみるとたいへん強力であった。例えば、NMR の残余双極子相互作用値を蛋白質の構造に当てはめる場合、400 残基分のデータをフィットするのに今の PC であれば 0.1 秒とかからない。当時でも 3 秒あれば十分であった。そこで、いろいろな初期値を 10 個ぐらい適当に準備して、それらを全部試してみる。 30 秒ぐらいで 10 個のフィッティングがすべて終わるので、その中から最小の二乗偏差を示す結果だけを選ぶ。エラーバーはどうするのかということになるが、これはモンテカルロ・シミュレーションがよい。100 個ほど模擬データを作っておき、全てにアメーバを振りかける。バナナを食べ終わる 10 分後にはすでに終わっており、その結果をエクセルで統計処理すればエラーバー（標準偏差）が出てくる。

結局つい数年前まで、NMR の緩和、残余双極子相互作用、回転拡散の異方性、拡散係数など、ほとんど全ての処理にアメーバを這わせたが、まったく問題なかった。緩和分散曲線もかなりは問題なし。もちろん Levenberg–Marquardt 法を使えば 0.1 秒ほどで終わってしまうが、今の PC をもってすると、それが延びたところで 1 秒ぐらいなので、もはや大した問題ではなくなっている。今は Perl に書き直したコードも使っている。C の 3-5 倍ぐらいかかるが、上記の NMR データの最適化で 1 分以上かかるような処理はあまりない。

(https://github.com/tikegami-tak/cpmg_hsqc) の中の perl/TI_lib_perl/TI_montecarlo.pm の中の lsquare_simplex というサブルーチンがアメーバ法です。

Fortran 一色という友達がおり「何故？」と尋ねてみると、いろいろ試してみたけれど、やはりこれが最速とのこと。特に数学処理では LAPACK ライブラリなどを使うと（もちろん C でも使えます）そこには何十年もかけて多くの人々が努力してきた成果が埋め込まれており、量子力学の計算ではそれを使わないと計算が終わるまでに一生かかるとのこと。NumPy も中ではこの LAPACK や BLAS といったライブラリーが使われているそうで、まあ道理で速いはずです。そこで、Perl でもこのようなライブラリーをさっと使えるようにならないだろうか？そしたら、アメーバからやっと卒業できるかも？一応、CPAN には Perl と LAPACK などをつなぐラッパーがあるようですが、NumPy のような使い方ではないような気がする。

ちなみに、最近はずるをして Mathematica で Levenberg–Marquardt 法をこっそりと使ったりするのであるが（Perl のアメーバには内緒）、Python を見てみると、変数の範囲を限定して Levenberg–Marquardt 法が使えるようになっていた。さらに CPMG 法では各残基ごとに異なる変数と分子全体で共通の変数があり、Mathematica のようにきっちりとまとまり過ぎた関数ではうまく最適化のためのコードを書けない時がある（ちゃんと勉強したら出来るのかもしれないが、私にはダメだった）。そんな変則的な時、ちょっと原始的な Python は非常に便利だった。また、アメーバ法も今はアルゴリズムがさらに進化しており、もっと高速化しているようです。

つれづれ思うこと Perl 編 その４

Perl と Python の比較だけでは話が単調になってしまうので、他の言語についても触れてみたい。

実は、最初に学んだ言語は Prolog だった。学んだというよりかは、大学の教養課程の授業「情報処理」で半年ほど触れただけである。1987 年当時、理系の授業で例として習うプログラム言語は普通は Pascal であった。ところが何故か私の授業の担当教員は Prolog に興味をもっており、そのため初めて触れるプログラム言語が Prolog となってしまった。「A さんの親は B さんで、B さんの兄は C さんで...」などという関係性を一行ずつ書いていく。そして、最後に A さんの叔父さんは誰ですか？」と質問すると「C さんです」という答が返ってくる。当時は人工知能のための言語として将来を期待されていたらしく、その先生も夢をもって語っていたのを思い出す。なにしろ初めて見たプログラムなので「Prolog で行列の計算ができないではないか！」などという疑問をもつこともなく、素直に「プログラミングって面白いな」と思った。

人口知能用の言語ということで、今の ChatGPT などの祖先のように思ってしまうが、おそらく全く違うように思う。今この言語はどうなっているのだろう？と思い調べてみたら、今もマニアの間では人気があるようだ。SWI-Prolog など無料でダウンロードできるサイトも活躍している。たぶん Prolog でしか処理できないような論理的な問題があるためであろう。C, Python, Perl などとは全く異なったコンセプトで設計されているので、視野を広げる意味でもいつか復習してみたい。なお、Perl との関係はなさそうで、唯一見つけられたのは、どちらもプログラムファイルの拡張子が ".pl" であることぐらいか？

Prolog の次に習ったのは C 言語で、これはたいへん難しかった。特にポインタの部分はよく分からない状態が続いた。しかし、当時勤めていた会社がたいへん教育的で、毎週 C 言語の講座を開いてくれて、少しずつではあるが理解できるようになってきた。1991 年当時は開発を担当していた半導体検査用電子顕微鏡を制御するのにアセンブリ言語（おそらくモトローラの MC6800 だったか?）が使われており、多くの技師がこれでプログラミングしていた。これは素人にはただの暗号に見える。そして、電顕の画像処理などに C 言語が使われ始めたのをきっかけに、機械制御の部分もワークステーション上での C 言語処理に少しずつ置き換えられつつあった。そのような状況もあり、アセンブリ言語の講座も社内にあった。しかし、その講座は一週間の合宿制になっており、参加するには上司の許可が必要であった。おそるおそる上司に相談してみると OK とのこと。実はこれが良かった。この講座によりアセンブリ言語がどのようにメモリーに値を入れたり、ビットをシフトさせたり、読み出したりしているかが分かり、同時に C 言語のアドレス（ポインタ）が理解できるようになった。あの時の上司の理解がなければ、プログラミングをその後 30 年以上にわたって今日までの仕事に活かすことはできなかったであろう。今も当時の教科書は残してあり、時々思い出しては感謝している（「マイコン回路の手ほどき」（1983）白土義男著）。

今はアセンブリ言語に触れることはあまりなく（マイコン制御ではあるのかな？）Python などあまりアドレスを意識しないで書けるプログラミング言語が主流となった。しかし、その根底にはメモリーの制御が常にある。よって、若い人ほどアセンブリ言語をちょっと勉強してほしい。私が受けた第二種情報処理技術者試験では CASL というアセンブリ言語が出題されたのであるが、今でもそれが選択肢のひとつとして残っていることを知って嬉しかった。是非 CASL を選んで受験して欲しいと思う。そしたら、人によっては、その後に学ぶすべてのプログラミング言語での見方が変わるだろう。

最後にアセンブリ言語と Perl との間に何か関係があるのか？という点になるが、もちろんあります！Perl でも Python でもそうですが、配列 abc をコピーしたりする時に (hiq = abc)、配列そのものが複製されるのか、それとも配列 abc を表す代表アドレスだけが別の変数 hiq に代入されるのかという問題を理解しておかないといけない。この仕様が言語によって異なるので大変ややこしい。さらに３次元 NMR スペクトルのように多次元データになってくると、配列そのものの複製なのか、それともアドレスだけの代入なのかの区別が複雑極まりない。C 言語では配列名がその配列の先頭アドレスを示すことになっているので、それをしっかりと理解しながらプログラミングできる。しかし、Perl にしろ Python にしろ他の言語になってくると、またその仕様が変わってくるので、いつもウェブで調べてはあれこれ悩んでいる（特に二次元目）。間違えて１万行 x １万行の行列を、何か要素を少しだけ計算するたびに何千回も複製していたら、メモリーがいくらあっても足りないですからね。。

まあ、しかし今の言語の良いところは、C 言語での malloc と free が不要な点でしょう。いわゆるメモリーの動的確保と解放です。この free を忘れると、演算をするたびに自由に使えるメモリー領域が減っていき、最後には PC が固まってしまう。このガーベッジコレクション（というのかな？）が裏で働いてくれるお陰で、安心してプログラムが作れるようになった。また、私が C 言語ではなく Perl に移った主な原因もこれです。

ちなみに、上記の電顕についていた画像処理装置（Cognex）は C 言語で動くのだが、フリーのメモリーの現在量を示してくれる関数があった。これは大変便利。自分の作ったプログラムを動かす前と後とでわずかでもメモリー量が異なった場合、それはどこかで free を忘れているということを意味する。もし、この機能がなかったら半導体工場のオートメーション上で数時間使ったら必ず再起動してしまうといった事故が起きていたことだろう。

ところで、うちは生物系であるので、学生達はあまりコンピュータについて習ってきていない。しかし、今後はスマホにしろ AI にしろ、コンピュータ技術と無縁ではいられないので、少しは知っておく必要がある。そこで、PC のメモリーを 8 GB から 16 GB に増設するついでに、PC の中身を開いて「これはハードディスク、これがシムメモリー、これが CPU」と、一応の基礎を教えることにした。そしたら、NMR の三次元スペクトルのシリーズを 10 個ぐらい開いてピークを拾い始めると、なぜ自分の PC が急に重くなったり、時には落ちてしまうのかの理由も分かるはず。ついでに、それらがスマホとどのように対応しているのかを知ると、ちょっと親近感が湧くようだ。そう、何でもブラックボックスとプロトコールだけで済ませず、その原理を知ると面白いものなのだ。

つれづれ思うこと Perl 編 その３

（その２）の続きとして、次は NumPy について。

例えば、筆者の NMR データの解析のためのプログラムには、exp(i H t) という計算が一杯入ってくる。ここで i は虚数、t は時間であるが、H はハミルトニアンと呼ばれる演算子であり、中身は行列である。つまり、e の何とか乗での冪数（べきすう）にあたる箇所に行列が来る。その行列も１万行 x １万列などという大きさの場合がある（NOESY, saturation transfer のシミュレーションなど）。3 の 2 乗は 9 とすぐに計算できるが、3 の行列乗（しかも内部には虚数）なんて、どのようにして計算するのだろう？と思われるかもしれないが、実はこの行列 (i H t) を対角化さえできれば簡単に計算できる。この exp(i H t) は物理では至るところに出てくるので、物理の計算では行列の対角化（つまり、固有値と固有関数の導出）ができるかどうかが重要な鍵であると言える。

これを Python で計算させると、リターンキーを押し込む前に速攻で答が返ってくる。NumPy と SciPy ライブラリーのおかげである。ここが Python の強いところ。このライブラリー部分は C 言語か Fortran かで書かれているらしいが、表面的にはまるで Python そのものが計算しているように見える。

一方、Perl でこれができるだろうか？一応 PDL（Perl Data Language）と呼ばれるライブラリーがあり、最近ではいろいろな計算ができるらしい。しかし、25 年ほど前は exp(i H t) が Perl で簡単に、しかもバナナジュースを飲み終える短時間内に計算できるかどうかはよく分からなかった。ここが個人的には勝敗を分けたような気がする。1996 年頃、何とかこの計算を Perl で行いたく、力づくでテーラー展開を使った。そうです、対角化という高等な操作を使わなくても、ひたすらテーラー展開をして汗水流せば、exp(i H t) を求めることができるのです。高校でよく x が小さい時 (1+x)^n = 1+nx ですよと習う、あの近似式もテーラー展開から来ている。これを１次でストップせずに、２次、３次と続けていくとますます真の値に近づいていくはずであるが、いったいどこまで続ければよいのだろう？

普通はまあ 10 次ぐらいまで計算して足し込めば十分だろうと思いがちである。ところが、どんどん次数を上げていくと、収束したり発散したりを繰り返したのである（行列の中に虚数が含まれているためだろうか？）。よって、変な次数のところでうっかり止めてしまうと、真値からかなり外れた値に落ちてしまう。というわけで、200 次ぐらいまで計算することにした。この 200 がよいのかどうかは今だに分からない。

ややこしそうに見えるが、要は行列の内積（行と列の要素を一つずつ掛け合わせてから足す）をどんどん繰り返していくだけである。これを計算するには、for 文による繰り返しをひたすら書いていけばよい。そして、NMR のシミュレーションでは FID 部分は 1,000 ポイントぐらいあるので、またこれを 1,000 回ぐらい計算するのである。いずれも for 文を使って。ちょっとややこしいのは、実数と虚数を分けたり混ぜたりしながら計算しないといけない点である（虚数 x 虚数 = 実数, 虚数 x 実数 = 虚数となる具合に）。

このようなプログラムを書いていると何故か汗が出てくる。計算するのはあくまで PC であって、自分が計算するわけではないのであるが、for 文をひとつ追加するごとに、まるで自分が激しく運動しているかのような錯覚に陥り、大変だなあと困って汗ばんでしまうのである。そして何とか書き上げ、いざ実行！だが、答が返ってこない。。。。バグっているのか、それとも PC がまだひたすら計算しているのかは？である。試しに１ポイントだけ計算させてみると、10 分ほど経って一応は正しい答が返ってきた。やった～と叫んだが、たった１ポイントだけではスペクトルをシミュレートするには程遠く、ほとんど役立たずである。それでそのまま放り出してやめてしまった。

月日は流れ 22 年ほど経ったある日、そのプログラムをふと思い出し、こわごわ動かしてみた。普通はそんな昔のプログラムが素直に動くことはない。しかし、この Perl はなにやら計算をし始めた。この Perl の後方互換性は大変すばらしい！私がこれまでに作った Perl の中でうまく動かなったのは一度だけで、それも 30 年以上前に作った初期版だけである。たしか、先頭にスペースのある行を split した時に、配列インデックスが 1 から始まったような気がする（今は０）。それで iMac-27 (core-i5, 32 GB memory）のマシンで 30 分後にプロンプトが返ってきた！恐る恐る gnuplot でグラフを表示させてみると、事前に Mathematica で計算したのと全く同じグラフが出てきたのである。これには大変驚きで、プログラムを前に「なんて可愛いヤツなんだ」と感動してしまった。

あまりもの嬉しさに論文に「Mathematica と Perl で図８をシミュレートした」などと書いたら、ある審査員 reviewer が「こいつはバカか？」と思ったのか、論文を reject すると同時に、NumPy と SciPy のプログラムを送ってきた。それでこれを実行してみると、3 秒ぐらいで同じグラフが画面に出てきたので（gnuplot もなしに）、共著者とともにまたまた驚き、なんて速いプログラムなんだ、大蛇、恐るべし！と、これまた感動してしまった。

それからは、もう勝てないものには勝てないと諦め、素直に Python で、しかし泣く泣く計算プログラムを作るようになった次第である。とは言え、もし Perl にもそのような高速の数学ライブラリが揃っていれば、さらにグラフも簡単に作れれば、間違いなく、蛇ではなくリャマをお友達にしている。

つれづれ思うこと Perl 編 その２

ソフトウェア言語を比べながら書くのであれば、もうちょっと段落ごとに内容をまとめるべきであるが、どうもそれは苦手であるので、頭に思い浮かんだことをそのまますぐにタイプすることにする。

そもそも何故、どのようなきっかけで Perl を始めたのだろう？

1994 年頃のことであるが、やはりデータ解析でファイルフォーマットの変換が必要になり、これは何かプログラムを使わなければと思った。しかし、その時に知っていたのは C 言語だけで、最初はプログラムをコンパイルして、出来上がった a.out を実行していた。しかし、PDB ファイルでアミノ酸残基番号を１だけずらすのに、C 言語でプログラムを書きコンパイルするなんて、ちょっと馬鹿げている（C ユーザの方、ゴメンなさい）。ある時、指導教員の机上に、微笑んでいるラクダの絵が描かれた分厚い本が置いてあり、この表紙がすこぶる気に入ってしまった。これが Perl を始めることになった唯一の動機かつ、きっかけである。その後、このシリーズの本には、ラクダ以外にも、ラマ、羊、山羊、豹などがあったような気がする。いま後ろの本箱を調べてみたら、リャマ、オオツノヒツジ、黒豹だった。リャマの絵などは大変かわいい。

一方 Python 本はというと、気持ち悪い大蛇！しかも、どこかで本物の蛇の写真が表紙に使われているのを見たことがある。これはとてもではないが、買って枕元に置こうとは思わない。ましてや学生に「この本で勉強してみたら？」と本を渡したとしたら、驚きようによっては嫌がらせかと誤解されかねない。ちなみに「Python １年生」という本が出ており、山羊、羊、犬などが登場していた。もちろん買いました！また、GAS プログラミングの本で「カエル」が表紙を飾っているのがありましたが、それもとても可愛かったので、学生が気に入っているようです。

この頃 1995、研究室の本棚には Python, Tcl/Tk, C++ の書籍もあり、今後どれを使っていこうか、どれを中心に勉強していこうかと迷った。その時に、Perl ではなく Python を選んでいたら、今頃ここにどのような文章を書いていたのだろうか？（C++ 本のリスも良かったのであるが、ちょっと怖い目をしていた。）

うろ覚えであるが、その頃から「Python は数学処理に向いている」とどこかで読んでそう思っていた。これが NumPy の始まりだったのだろうか？（初版は 1995 年らしい）。指導教員とともに一か月ほど使ったことがあるが「やはり Perl の方が良いよね」ということで、また Perl に戻ってしまった。その頃は Python の影は薄く、Perl とそっくりな言語がまた作られたぐらいにしか思っていなかった。

その後もずっと Perl の方が人気が高かったように思うが、いつ頃からか 2010 年頃？AI の台頭とともに急に Python が登ってきた。今ではなんでもかんでも Python になって、某雑誌なんて何か月間も Python 特集が続いている。この優劣を分けた鍵は何だったのだろう？

私が個人的に思う鍵は、NumPy の存在である。これについては、次の（その３）で書いてみたい。それに加えて、ソフトウェア言語の人気には、予測の難しい一種の非線形的な現象が絡んでいるように思う。ある言語の人気が高まれば、その言語でいろいろなライブラリーが拡充してくる。すると、ますます多くの人が使うようになる。逆にその他の言語はますます使われなくなる。こうして、あっという間に圧倒的な差が生じる。このような現象は協同的といえる。

生物の世界でもこの協同現象は至るところに見つかる。例えば血中のヘモグロビンなどは、４つのサブユニットから構成されているが、一つ目のサブユニットに酸素がつくと、次のサブユニットには酸素がより付きやすくなる。すると、三番目にはますます酸素が付きやすくなる。逆に酸素が離れ始めると、ますます他のサブユニットからも酸素が外れやすくなる。この「more and more」「less and less」の現象により、肺では４つのサブユニット全てに酸素がくっつき、一方、体の端の毛細血管の中では酸素を全て組織に渡すことになる。このようにヘモグロビンの四量体構造はスィッチのような機能を発揮するのである。その他にもプロウィルスの活性化など、生物における協同性には枚挙にいとまがない。

いったんシーソーが傾き始めると、雪崩のように、もう止めることができなくなる程に独りでに傾いていく。よって、Python の人気にもこの協同性が生じてしまったのではないかと思う。その最初の一押しが NumPy だったのではないだろうか？

協同性には怖い裏面もある。先ほどのヘモグロビンのように逆行もまた可能だからである。つまり、いったん人気を失い始めると、まるでスィッチを切るかのように急降下する。しかし、Julia なども思ったほど伸びてきていない現状で、今のところまだ大蛇王国が続きそうな気がする。

さて、急降下してしまった Perl はどうなるのだろう？Perl7 が出てくれればと願っていたのに、なんだかまた降下してしまった。もう別に Perl5 のコピーそのままでよいので、とりあえずは Perl7 という旗を揚げてしまえばよいのではないだろうか？まあ、人気は落ちるだろうけど、きっと誰かが管理しながら、なんとかいつでも使える状況があと何十年かは続くことを期待している。

つれづれ思うこと Perl 編 その１

あまり堅苦しい論文の話ばかりを書いていても一向にブログが進まないので、もうちょっと気軽な内容も書くことにしました。

蛋白質の立体構造の計算でいつも困るのは、Protein-Data-Bank の昔のフォーマットが統一されていないこと。例えばアミド水素など、あるソフトウェアでは 'H' と書かれていないといけないのに（おそらくこれが正式）、別のソフトでは 'HN' と書かれていないと動かない。その他にも CH2 の β1, β2, β3 の割り当てなど、泣きたくなるような不一致がいくつもある。これを修正するのは一苦労で某友人もイライラが募るとのこと。なお、今後は mmCIF と呼ばれるフォーマットに統一されていくので、このような問題はなくなるようです。

そういう時にファイル変換のためのプログラムが役立つわけであるが、どの言語を使おうか？Python, Perl, Ruby, JavaScript どれでもよいのであるが、そもそもプログラミングが本職ではないので、できるだけ手っ取り早く、さっさとファイル変換を済ませたい。すると、どうしても過去に自分で作ったプログラムを少しだけ修正して目前の問題に間に合わせることになる。これが 30 年間も続くと、ずっとその言語一色という状況が起きてしまう。

それで私は Perl 一色だったのであるが、最近はウェブで「Perl 将来性」などと引くと、「オワコン = 終わったコンテンツ」などと嘆かわしいページばかりが目につくようになった。それとは対照的に Python はまだまだ鰻上りである。というわけで、数年前からちょっとでも慣れようと思い "無理して" Python で書くようにした。Perl も Python もどちらも似たものだろうと舐めてかかったのであるが、なんと Perl の何十倍も作るのに時間がかかってしまった。

何故か？決定的だったのは括弧の有無である。例えば Perl の場合

if ($abc == 0)
{
     $def = 3;
}

のように {  } を使うが、Python では使わずに行を右へずらす（インデント）。ちなみに、私は癖で

if ($abc == 3) {

    $def = 3;
}

のようには書かない。その理由は、{ と } がエディター上で縦に揃っていると、if 文の及ぶ範囲が一目で分かるためである。しかし、Python では括弧がないので、vi で Python プログラムを編集した時、if 文の階段が５段ぐらい入れ子になり、そこでインデントを間違えた瞬間に崩壊してしまった。間違えた！と気づいても、どこをどう修正してよいのか分からなくなるのである。他人の Python プログラムを破壊するのは簡単である。エディター上でスペースキーをちょこっと押してやると、もう動かなくなり、修正にたいへん時間がかかる。

一応、vscode などの高機能エディターを使うと、縦の補助線が出るので少しは助かる。しかし、昔の人はこのような補助線なしにどうやって Python プログラムを書いていたのだろう？

Perl にも面白いところがあって、例えば "123" などを数値とみてもよいし、文字列と見てもよい。計算の中でその変数を使うと（例えば、1 を足すなど）自動的に数値として解釈してくれる。しかし、Python ではそうではないので、ちゃんと（昔習った C 言語のように）数値と文字列の違いをきっちりと区別し、ある時にはいずれかへの変換をちゃんと書いてあげないといけない。「そんな事は当たり前だろう」と言われそうであるが、30 年間も Perl 温泉に浸かってしまうと、なんて Python は不便なんだと悪口を吐いてしまう。。。

あまり Perl ばかりを褒めていてもよくないので、Python を使い始めて「これは良い」と思った点も書こう。それは、Perl では変数はデフォルトで大局変数（グローバル）になっているのに対して、Python ではデフォルトで局所変数（ローカル）になっている点である。大局変数の方が書きやすいのであるが、時に変数名が重複してしまい、知らない間に保存してあった変数の値が変わってしまうこともあり得る。そのため、Perl では、これはローカル変数にしたいと思った時に変数の前に my を付ける。一方、Python では考え方が真逆であるので、普通に変数を設ければ、それは自動的に局所変数となる。よって、安心してプログラムを書ける。しかし、これも 30 年間グローバル変数に慣れてしまうと、Python で書いた時に「なんて窮屈なんだ。この何十個もある変数をどうやって子供のサブルーチン関数に渡したらよいのか？まさか全部を引数で？」と悩んでしまう（もちろん、オブジェクトにすればよいのだけれど）。

これを書いていて、今もう一つ Python で良かった点を思い出した。それは、引数にデフォルト値を設定できることである。これは便利。もしかして Perl にもちゃんと備わっていて、私が知らないだけなのだろうか？