雑感等

音楽,数学,語学,その他に関するメモを記す.

Vossのアルゴリズムとシェンカー分析

Five Graphic Music Analysesの"Ich bin's, ich solle büssen"の分析を見た.

Url.Tafelから3.Schicht, 2.Schicht, 1,Schicht, Ursatzと見比べていくと,
VossのアルゴリズムでfBm信号を生成する過程の逆に似ていると思った.

VossによるfBm生成過程は,あるサンプルと直後のサンプルとの中点に乱数を加算し,サンプルを増やしていくことだ.
シェンカー分析において表層から深層にたどる過程は,増えていったサンプル(=労作された音符)を取り除いていくことだ.

旋律生成のために頭音を労作する方法が同書の分析から読み取れる.

  • 頭音(3^, 2^, 1^)を二点に配置し,その二点を滑らかにつなぐ.
    • 例1:"Ich bin's, ich solle büssen"のような場合
    • e(3^) e(3^):二点に配置する.
    • →e(3^)-d-c-d-e(3^):滑らかにつなぐ.この場合はe-(三度下降)-c-(三度上昇)-e
    • 例2:BWV 846のような場合
    • e(3^) e,(3^):二点に配置する.二つ目のeはオクターブ下
    • →e(3^)-d-c-h-a-g-f-e(3^)

以上のような労作の方法は『調性音楽のシェンカー分析』に網羅的に示されていたが,感覚的にわかりにくかった.
Five Graphic Music Analysesでは,1曲全体を複数の構造レヴェルで分析した結果が示されており,感覚的にわかりやすい.

文献

  • H. Schenker and F. Salzer: Five Graphic Music Analyses, 2nd ed. (Dover Publications Inc., New York, 2012).
  • アレン・キャドウォーラダー, デイヴィッド・ガニェ,角倉一郎訳: 調性音楽のシェンカー分析. (音楽之友社, 2013).

文献の示し方

出典:若林敦. 理工系の日本語作文トレーニング. 朝倉書店, 2000.

(i) 図書を示す場合([ ] は書かない場合もある):

著(編)者名 書名 [版表示] 出版社(者)またはシリーズ記載 [出版地] 出版年 [ページ]

[版表示]=初版は不要.改訂版・増補版及び第2版以後の版は必ず表示する.
[ページ]=引用箇所,参照箇所を特に指定する場合に表示する.
[出版地]=日本語の文章では,書かないことが多い.

(ii) 図書収蔵論文を示す場合:

論文の著者名 論文名(表題) 所収図書の編(著)者名 書名 …途中は,上の図書を示す場合に同じ… [ページ]

[ページ]=記載する場合は所収論文の掲載されているはじめのページとおわりのページを記す.おわりのページは省略してもよい.

(iii) 雑誌論文,雑誌記事を示す場合:

著者名 記事タイトル,論文名(表題) 誌名 巻(号) 発行年 [ページ]

[ページ]=記載する場合は論文,記事の掲載されているはじめのページとおわりのページを記す.おわりのページは省略してもよい.

フラクショナルブラウン運動の生成するmatlabプログラム(Vossのアルゴリズムによる)

function sig=voss(hurst,len)
%Vossのアルゴリズムによってフラクショナルブラウン運動を生成する.
%husrt:Hurst指数
%len:出力するデータ長

%2^n+1>=lenとなる最小のnを求め,その長さでfBmを生成する
lenFbm=2^ceil(log(len-1)/log(2))+1;

%初期化
fBm=zeros(lenFbm,1);

%fBmの終点と始点の値を決める.
fBm(lenFbm)=0;
s2_0=1; %fBm(lenFbm)の値を決める乱数の分散.このパラメタによりfBm全体の振幅が変わる.
fBm(lenFbm)=randn()*sqrt(s2_0);

%再帰的に残りの点を決める.
fBm=calcmidpoint(1,1,lenFbm,fBm,hurst,s2_0);

%出力するデータの長さを指定された長さで打ち切る
sig=fBm(1:len);

end

function fBm=calcmidpoint(n,a,b,fBm,hurst,s2_0)
%n:再帰の深さ
%a:区間の左端要素の添字
%b:区間の右端要素の添字
%fBm:フラクショナルブラウン運動の時系列データ
%hurst:ハースト指数
%s2_0:fBm(end)の値を決めるのに使った乱数の分散

%中点の添字を計算
mid=(b-a)/2+a;

%中点が整数番目の添字を持つなら中点の値を計算する.
if mid==fix(mid)
    
    %s2:分散
    s2=( ( 1-2^(2*hurst-2) )*s2_0 )/( 2^(2*hurst*n) );
    
    %中点の値=区間の両端の値の平均値+分散
    fBm(mid)=(fBm(b)+fBm(a))/2+randn()*sqrt(s2);
    fBm=calcmidpoint(n+1,a,mid,fBm,hurst,s2_0);
    fBm=calcmidpoint(n+1,mid,b,fBm,hurst,s2_0);
end
end

Vossのアルゴリズムに関する式変形の過程(フラクショナルブラウン運動の生成)

文献: 本田勝也. フラクタル. 朝倉出版, 2002. pp. 110-112.

ラクショナルブラウン運動の関数を X(x),\, (0\leq x \leq1)を生成する「ヴォスのアルゴリズム」が文献に示されていた.
その中の式変形を詳細に示す.

 X\left(  0 \right) =0,\, X\left(  1 \right) =\delta_{0}, \, (\delta_{0} \sim \text{平均0,分散1の正規分布})
 X \left( \frac{1}{2} \right) = \frac{1}{2} \left\{ X(0)+X(1) \right\} + \delta_{1},\,\left( \delta_{1} \sim \text{平均0,分散$\sigma_{1}^{2}$の正規分布} \right) \hspace{3em}\text{(10.17)}
ただし X(1)-X(0) \delta_1は互いに独立である.

 \langle \left\{  X \left( \frac{1}{2} \right) - X \left(  0 \right)   \right\}^{2} \rangle \hspace{3em}\text{(10.18)左辺}
式(10.17)を代入
 =\langle \left\{  \frac{1}{2} \left\{ X(0)+X(1) \right\} + \delta_{1} - X \left(  0 \right)   \right\}^{2} \rangle
展開
 =\langle \left\{  \frac{1}{2}  X(0) + \frac{1}{2} X(1)  + \delta_{1} - X \left(  0 \right)   \right\}^{2} \rangle
項をまとめる
 =\langle \left\{   \frac{1}{2} X(1) -\frac{1}{2}  X(0)  + \delta_{1}    \right\}^{2} \rangle
 =\langle \left\{   \frac{1}{2} \left\{ X(1) - X(0) \right\} + \delta_{1}    \right\}^{2} \rangle
二乗を展開
 =\langle   \frac{1}{4} \left\{ X(1) - X(0) \right\}^{2} + {\delta_{1}}^{2} +  {\delta_{1}} \left\{ X(1) - X(0) \right\} \rangle
加法性を適用(統計平均と統計の期待値(平均)は同じ?)
 =\langle  \frac{1}{4} \left\{ X(1) - X(0) \right\}^{2} \rangle + \langle {\delta_{1}}^{2} \rangle + \langle  {\delta_{1}} \left\{ X(1) - X(0) \right\} \rangle
 X(1)-X(0) \delta_1は互いに独立だから,積の平均を平均の積に分解
 =\langle  \frac{1}{4} \left\{ X(1) - X(0) \right\}^{2} \rangle + \langle {\delta_{1}}^{2} \rangle + \langle  {\delta_{1}}\rangle  \langle  \ X(1) - X(0) \rangle
正規分布に従う乱数の平均だから   \langle {\delta_{1}}^{2} \rangle = \sigma_{1}^{2}, \,   \langle  {\delta_{1}} \rangle= 0
 =\langle  \frac{1}{4} \left\{ X(1) - X(0) \right\}^{2} \rangle + \sigma_{1}^{2} + 0 \times  \langle  X(1) - X(0) \rangle
 =\langle  \frac{1}{4} \left\{ X(1) - X(0) \right\}^{2} \rangle + \sigma_{1}^{2} \hspace{3em}\text{(10.18)右辺}


他の書籍を読んだ際,統計平均を\langle  \rangleで表したり,ほぼ等しいことを \sim で表したりする記法を知らなかったため,式が理解できないことがあった.
専門分野の表記はよそ者がみると理解できない.

フラクショナルブラウン運動の分散特性のプロット

matlabプログラムで分散特性\sigma^{2}_{H}= \langle \lvert f(t+\tau)-f(t) \rvert ^{2} \rangleに対する時間間隔 \tauを両対数グラフにプロットする.

Hurst指数がHのフラクショナルブラウン運動(fBm)の信号に対して分散 \sigma^{2}_{H}-時間間隔 \tauをプロットすると,
両対数グラフ上で傾きが2Hの線形となる範囲が現れる.

フラクショナルブラウン運動の生成プログラム - 雑感等に記載したプログラムを用い,
Hurst指数H=0.2, 0.5, 0.8のfBm信号を生成し,分散 \sigma^{2}_{H}-時間間隔 \tauをプロットした.

 \tauが小さい領域では線形にプロットされる.

H=0.2
f:id:kazmus:20180804111504p:plain

H=0.5
f:id:kazmus:20180804111517p:plain

H=0.8
f:id:kazmus:20180804111529p:plain

分散 \sigma^{2}_{H}-時間間隔 \tauをプロットするmatlabプログラム

%使用例
%fBm=csvread('fBm.csv');
%dif_time(fBm)
%第二引数は空で良い

function dif_time(sig,fs)
if 1~=exist('fs')
    fs=length(sig);
end
ts=1/fs;
loopMax=length(sig);
meanMoment=zeros(loopMax,1);
subplot(211)
for k=1:100:loopMax
    diffSquareSig=(sig(k+1:end)-sig(1:end-k)).^2;
    loglog((k*ts), ( mean(diffSquareSig)) ,'o')
    hold on
end
ylabel('<|f (t+\tau)-f (t)|^{2}>','FontSize',16)
xlabel('\tau','FontSize',16)
subplot(212)
plot(sig)
ylabel('amplitude','FontSize',16)
xlabel('time','FontSize',16)
set(gca,'FontSize',12)
end

フラクショナルブラウン運動の生成プログラム

正規分布に従う乱数を分数階積分してフラクショナルブラウン運動を生成し,csv形式に出力するプログラム.
下記のC++11のコンパイラで動作を確認

  • Embarcadero Technologies Inc. bcc32c version 3.3.1
  • Microsoft Visual Studio Community 2015 Version 14.0.25431.01 Update 3
//>>bcc32c fBm_generator.cpp
//>>fBm_generator.exe 
//または
//>>fBm_generator.exe [Hurst指数]
#include<iostream>
#include<math.h>
#include <random>
#include <fstream>
#include <stdlib.h>

////////設定項目////////

//Hurst指数 下記の値または,コマンドライン引数の一番目
#define HURST_EXPONENT (0.5)
//出力する時系列データの点数
#define DATALEN (65536)

//fBmの時系列データを出力するファイル名.ただし,出力形式はcsv
#define fBm_FILENAME ("fBm.csv")

//fBmの生成に使用する乱数の時系列データを出力するファイル名.ただし,出力形式はcsv
#define WHITE_NOISE_FILENAME ("white.csv")


////////計算プログラム////////

//分数階積分する関数"FractionalIntegral"は下記のサイトから引用,一部改変した.
//http://www.geocities.co.jp/Technopolis-Mars/1795/study/fractint.html
#define MAXKAISUU s //積分演算をおこなう時に利用する過去の情報の最大数
#define m_chk_mean false //最終的に平均値を零とする処理をおこなうかどうか
#define m_DT 1.0 //時系列のサンプリング時間間隔

//入力信号をそのまま微分する
//*data・・・微分をおこなう対象データのポインタ
//sekibun・・・積分階数
//s・・・データの大きさ
void FractionalIntegral(double *data, double sekibun, int s)
{
	register int m, t;

	//sekibunが整数かどうか判別する(整数の時には動作が異なる)
	int kaisuu = (sekibun == (int)sekibun) ? (int)sekibun : MAXKAISUU;
	double* ganma = new double[kaisuu + 1];

	//ガンマ関数の値をストックする    (n+1)C(m+1)
	//ついでに符号も保存する
	ganma[0] = 1.0;
	for (m = 1; m <= kaisuu; m++)
		ganma[m] = -ganma[m - 1] * (sekibun - m + 1) / (double)m;

	//非整数階積分をする
	double p = pow(m_DT, sekibun); //h^d
								   // for(t = s; t; data[--t] -= data[0]);  //積分の場合はいらない

	register double f;
	int M;
	for (t = 0; t < s; t++) {
		M = (t < kaisuu) ? t : kaisuu;
		for (f = data[t] * p, m = 1; m <= M; f -= ganma[m] * data[t - m], m++);
		data[t] = f;
	}
	delete[] ganma;

	//データの平均を0に
	if (m_chk_mean) {
		for (f = 0.0, t = s; t; f += data[--t]);
		for (f /= s, t = s; t; data[--t] -= f);
	}
}

//正規分布に従う乱数を発生させる関数"randnormal"は下記のサイトを引用,改変した.
//https://cpprefjp.github.io/reference/random.html
double randnormal(double mean, double stdDev) {
	// メルセンヌ・ツイスター法による擬似乱数生成器を、
	// ハードウェア乱数をシードにして初期化
	std::random_device seed_gen;
	std::mt19937 engine(seed_gen());

	// 正規分布
	// 平均mean、標準偏差stdDevで分布させる
	std::normal_distribution<> dist(mean, stdDev);
	return dist(engine);
}


int main(int argc, char *argv[]) {
	double signal[DATALEN];

	//正規分布に従う乱数を発生
	for (int i = 0; i < DATALEN; ++i) {
		signal[i] = randnormal(0, 1);
	}
	std::cout << "乱数の一番目の値:" << signal[0] << std::endl;
	//正規分布に従う乱数を保存
	std::ofstream filew(WHITE_NOISE_FILENAME);
	for (int i = 0; i < DATALEN; ++i) {
		filew << signal[i] << "\n";
	}

	//Hurst指数に0.5を加算し,乱数時系列信号を分数階積分する階数integrateOrderを決定
	double integrateOrder = HURST_EXPONENT + 0.5;

	if (argc > 1) {
		char *strtod_error;
		double input = strtod(argv[1], &strtod_error);
		if (*strtod_error != '\0') {
			std::cout << "Hurst指数として設定できない文字列:" << strtod_error << std::endl;
		}
		integrateOrder = input + 0.5;
	}
	std::cout << "Hurst指数:" << integrateOrder - 0.5 << std::endl;
	//分数階積分を実行
	FractionalIntegral(signal, integrateOrder, DATALEN);

	//フラクショナルブラウン運動信号を保存
	std::ofstream file(fBm_FILENAME);
	for (int i = 0; i < DATALEN; ++i) {
		file << signal[i] << "\n";
	}

	file.close();
	filew.close();
	return 0;
}