LISP

信号の包絡線を計算する (Common Lisp)

信号の包絡線 (signal envelope) は「信号の頂点にシーツをかぶせたときのような概形」を持つ波形です。 MatlabやRで書いた信号包絡線の計算プログラムをCommon Lispに移植しました。信号包絡線の計算方法は振幅の絶対値に対して、(1)低域フィルタをかける、…

ピンクノイズを生成するLISP関数

以前PythonやJSFX用のものも作ったのですが、ピンクノイズはいろいろなところで活躍するので、Common Lisp用にも作っておきます。 ピンクノイズを作るPython関数 - 丸井綜研 REAPERのJSFXでピンクノイズを生成してみた - 丸井綜研 Voss-McCartney法 DSP gene…

重回帰分析 (Lisp Advent Calendar 2023)

Lisp Advent Calendar 2023の第1日目の記事です。Lispについては全くのヘナチョコなので1日目は恐縮なのですが、カレンダーの先頭が空になっているのも幸先が悪いので……。 さて、車輪の再発明を続けて、ようやく逆行列の計算までたどりつきました。今回は逆…

逆行列の計算

Common Lispで逆行列の計算をするための準備をしてきました。 行列式の計算(Common Lisp) - 丸井綜研 LU分解(Common Lisp) - 丸井綜研 LUP分解(Common Lisp) - 丸井綜研 連立方程式を解く(Common Lisp) - 丸井綜研 今日はついに逆行列の計算です。と…

連立方程式を解く

先日以下の記事にしたように、LUP分解をすると正則行列𝐴に対して𝑃𝐴 = 𝐿𝑈という分解ができるのでした。𝐿は下三角行列、𝑈は上三角行列、𝑃は並べ替えの情報が入っています。 marui.hatenablog.com これを使って連立方程式𝐴𝑥 = 𝑏を解くプログラムを作ります。LU…

LUP分解

【2023-11-09:入力行列を書き換えながら処理をするようになっていたので、呼び出し元でも行列の内容が書き変わってしまっていました。複製を作って作業するように修正しました。】 昨日LU分解のコードを掲載しましたが、単純なアルゴリズムなので、𝑎₁,₁要素…

LU分解

10月29日(日)にShibuya.lispのもくもく会がありました。残念ながら雑談の途中で中座しないといけなくなってしまいましたが、とても楽しかったです。この日のもくもく会の目標は「逆行列を計算できるようにすること」でした。実際は途中までしかできなかった…

行列式の計算

Common Lispで車輪の再発明をし続けています。いろいろと統計に使う関数を作っていましたが、そろそろ逆行列を計算したくなってきました。逆行列と固有値は線形代数の基本ですものね。*1 逆行列の計算をする前には、行列式(determinant)という、逆行列が計…

コラッツ問題

コラッツ問題という数学の未解決問題があります。自然数nが偶数のときn÷2を計算、nが奇数のとき3n+1を計算し、その結果に同じルールを適用し続けると、最終的に1に至る(あるいは最終的に1→4→2→1の無限ループになる)、というものです。すべての自然数につい…

cl-portaudioを使った音の入出力

前回の記事でCommon Lisp上でPortAudioが使えるようになりました。 marui.hatenablog.com 今回は音の入出力をやってみます。基本的にはcl-portaudioのドキュメントに載っている「入力と出力を15秒間だけつなげる」というデモプログラムを実行するだけです。…

LispからPortAudioでデバイス情報を得る

Common Lispからcl-portaudioを使ってPortAudioを呼び出して音を再生したいのですが、まずはデバイス情報を得るところまでやりました。 cl-portaudioを利用するには例によってQuicklispを使うと便利みたいです。 (ql:quickload :cl-portaudio) とはいえ、bod…

配列の環状シフトなど

昨日「相互相関の計算 - 丸井綜研」で「circshiftのあたりが面倒な気がしたので」と書きましたが、そんなに面倒ではなかったので書き残しておきます。 (defun circshift (arry k) "Circular shift the array ARRY by amount K. Shift towards right if K is …

相互相関の計算

相互相関 (cross-correlation) を計算するCommon Lispプログラムを作りました。とは言っても、畳み込み演算ができれば簡単に実装できます。Matlabのxcorr関数とそれなりの互換性があるものが欲しかったので、相互相関値だけでなくラグも返すようにしました。…

t検定(1群)

「二項検定と二種類の過誤 - 丸井綜研」で、Common Lispで二項検定ができるようになりました。今日は1群のt検定を実装したいと思います。 t分布 t検定をするには、まずt分布の確率密度関数と累積分布関数を作る必要があります(検定だけなら累積分布関数だけ…

ガンマ関数(再訪)

6月25日(日)にShibuya.lispのもくもく会がありました。14:30にオンラインで集合して、自己紹介と質問タイム(勉強になった!)があり、15時頃~17時頃にもくもくと作業をしていきました。(もくもく会の後の雑談は予定を1時間近くオーバーして盛り上がりまし…

applyとreduceの違い

リスト内の数値の合計を計算しようとして (defun sum (lst) (apply #'+ lst)) (sum '(1 2 3 4 5)) ;=> 15 というような関数を書きました。普通に使うぶんには全く問題なかったのですが、長いリスト(具体的には480,000個のdouble-float)を引数にしたら「ス…

二項検定と二種類の過誤

先日、二項分布の累積分布関数を計算できるようになりました。 marui.hatenablog.com これを使って、簡単な二項検定をやってみたいと思います(あまりLISPとは関係ない内容ですが)。題材としては、ジャンケンにやたら強い人がいて「俺には天から与えられた…

二項分布の累積分布関数(と正則化不完全ベータ関数)

統計に使える関数の実装を続けています。2月のLISPもくもく会では二項分布の確率質量関数(dbinom)を実装しました。今度は累積分布関数(pbinom)が欲しいところです。 marui.hatenablog.com 原理どおりに動くのであれば確率質量関数を累積加算すると累積分…

二項分布の確率質量関数(改)

以下のエントリに、二項分布の確率質量関数をCommon Lispで計算する関数を紹介したのですが、どうも意図した通りに動いてくれないときがありました。少し時間を作って調べてみたら、コードに間違いがあるのが明らかになりました。 marui.hatenablog.com Load…

二項分布の確率質量関数

【2023-04-17】このページの最下部に載っているコードでは、特定の引数の時に間違った答えが出力されたり、無限ループに入ってしまうことが分かりました。たとえば(dbinom 4 10 0.5) ;=>0.557となってしまいます(正しくは0.205)。元のアルゴリズムには問題…

PGMファイルへの画像書き出し

車輪の再発明を続けています。今日は画像の書き出しです。 スペクトログラムを計算できるようになったものの、それを画面表示するためにはなんらかの画像にしないといけないんですね。画面にウィンドウを作って表示するよりもファイルにするほうが簡単なので…

スペクトログラム

Common Lispで振幅スペクトルの計算ができるようになり、Hann窓もかけられるようになったので、次はスペクトログラムの実装にいきたいと思います。 スペクトログラムは、振幅あるいはパワーの強さを時間と周波数ごとに計算したものです。以下の画像ではAudac…

ガンマ関数

非負整数nの階乗n!はn×(n-1)×(n-2)×(n-3)×...×2×1と定義されています。また、0!=1とします。例えば6!であれば、6×5×4×3×2×1=720です。これを、整数以外にも拡張したのがガンマ関数で、nが整数の時にはΓ(n+1) = n!となります。例えばΓ(6)=120、Γ(7)=720となり…

正規分布の確率密度関数と累積分布関数の計算

目次 正規分布の確率密度関数 累積分布関数の近似式(Abramowitz and Stegun) 累積分布関数の逆関数の近似式(Abramowitz and Stegun) 累積分布関数と逆関数を誤差関数から計算する(Numerical Recipes) Common Lispを音だけでなく、統計にも使いたいです…

フィルタの可視化(freqz相当の関数を作る)

信号にフィルタをかけられるようになり、そして、レシピ通りのフィルタを作れるようになりました。 IIRフィルタをかける関数(Lisp Advent Calendar 2022) - 丸井綜研 双二次フィルタの係数(Lisp Advent Calendar 2022) - 丸井綜研 次に、フィルタの周波…

双二次フィルタの係数(Lisp Advent Calendar 2022)

この記事はLisp Advent Calendar 2022の25日目です。Lispで音プログラミングに向けた、さらなる一歩。 目次 はじめに RBJ Cookbookによる低域フィルタ フィルタを使ってみる 他の形状のフィルタ はじめに 先日フィルタをかける関数を作りました(IIRフィルタ…

IIRフィルタをかける関数(Lisp Advent Calendar 2022)

この記事はLisp Advent Calendar 2022の16日目です。Lispで音を自由自在に扱うという野望に向けての次の一歩。 デジタルフィルタ ちょっと前に高速フーリエ変換を用いて畳み込み演算の実装をしました。2信号の時間領域での畳み込み演算は周波数領域での要素…

一様乱数と標準正規乱数の行列を得る(Lisp Advent Calendar 2022)

この記事はLisp Advent Calendar 2022の13日目です。 目次 はじめに 一様乱数行列 改良 標準正規乱数(おまけ) 標準正規乱数行列について追記(2022-12-18) はじめに 任意の大きさの行列に乱数が入っているものが欲しいことがあります。Matlabのrand(5, 4)…

音ファイルを書き出すwavwrite関数(Lisp Advent Calendar 2022)

この記事はLisp Advent Calendar 2022の11日目です。Lispには入門と挫折を繰り返していましたが、今年4月にLispに再入門し、なんとか日常的なイロイロに使えるようにしているところです。 普段使いするためには自分の専門分野のプログラムを作るのがよいだろ…

Common Lisp用のwavread

(追記:関連記事として「音ファイルを書き出すwavwrite関数(Lisp Advent Calendar 2022) - 丸井綜研」を書きました) 9月に「Common Lispで音ファイルを読み込む - 丸井綜研」というエントリでcl-wavとbodge-sndfileを試用しました。Apple Siliconではbod…