非整数ブラウン運動の説明用のアニメを作ってみました． ブラウン運動 (H = 0.5)のサンプルパスの拡大非整数ブラウン運動 (H = 0.8)のサンプルパス非整数ブラウン運動 (H = 0.3)のサンプルパス高解像度な，Youtube版は以下です． www.youtube.com www.youtub…

Detrended Fluctuation Analysis (DFA)で，ホワイトノイズを解析したらどうなるか，数値実験で確かめてみます．理論的な結果は，前の記事を参照してください． chaos-kiyono.hatenablog.com 数値時系列の生成 1次DFA 2次DFA 数値時系列の生成 Rでホワイトノ…

Detrended Fluctuation Analysis (DFA)で，無相関な信号，つまり，ホワイトノイズを解析したらどうなるかのお話です．DFAに詳しい人の多くは，そんなのスケールをとして，ゆらぎ関数 (fluctuation function)は， だろというかもしれません．が大きい，漸近的…

昨日，名古屋に出張しました．新幹線は，コロナ前と比べればまだ空いていますが，海外の旅行者も戻ってきて，席が埋まってきたように感じます．ということで，今回は，Rで乱数を生成するコマンドの整理です．疑似乱数を使って下の図のような，中心極限定理の…

世の中には正規分布とか，指数分布とか，ワイブル分布とかいろいろな分布に従う現象があります．確率的な現象が従う分布の関数形を特定すると，背後にある生成メカニズムを考えるヒントになるので，分布のことを理解しておくと役に立ちます．今回は，Rで確率…

私は肉まんの白い部分が好きです．ほのかに甘い，あの白い部分だけ食べたいです． ということで，今回は最近のFFTのお話 chaos-kiyono.hatenablog.com の続きです． FFTの説明で，バタフライ演算ってのが登場します．バタフライ演算は下の図のようなバタフラ…

前回，2のべき乗の長さの時系列を，高速に離散フーリエ変換 (FFT)できるアルゴリズムのRスクリプトを作りました (高速じゃないけど)． chaos-kiyono.hatenablog.com このFFTアルゴリズムは，Cooley-Tukey型と呼ばれるそうです (知らんけど)．今回はこのアル…

Rとか，Pythonとかを使って，誰でも簡単に高速フーリエ変換 (fast Fourier transform, FFT)を使える時代になりました．最近コロナで，人と会わなくなって，久しぶりに遠隔会議であった方々が老けたなーと感じることが多くなりました．同様に私も年寄りになっ…

プログラムを作って何か処理するとき，処理時間が短い方が良いに決まってます．ということで，今回はRの処理時間についてのお話です． 処理時間を無駄に長くしないための一般論として，私が説明できることは，「forループは使わずにベクトルのまま計算した方…

1時間ほど前に書いた記事では，Rで文字列を時刻データに変換するのに「私はas.POSIXctを使います」と言ってましたが，訂正します． 私はas.POSIXctを２度と使いません． 今なら断言できます．Rでas.POSIXctを使ってる奴，素人です． as.POSIXctを使ってる人…

私はサツマイモ好きなので，コンビニのサツマイモ系スイーツは必ずチェックします．今日は大学のローソンで「さつまいもあんぱん コグマパン」という菓子パンが新発売されていたのでさっそく買って食べてみました．ということで，今回はRで時刻を日付に変え…

オムロン活動量計でHJA-750Cというのがあります．www.healthcare.omron.co.jpこれは研究や臨床検査用なので，一般には発売されていません．大学などの研究機関に所属していれば購入可能です． といことで，今回はHJA-750Cを使って研究をしたい人向けのに，R…

最近はYoutubeでマリマリマリーの動画を見ています．その影響で，ひ○ゆきのものマネができたらいいなーと思うようになりました．皆さんも，次のセミナーではひ○ゆき風のしゃべりでお願いします． ということで，今回は前に紹介したコサイナー解析 chaos-kiyo…

今回はファイルからデータテーブル (表)として読み込んで，データを処理する方法をいくつか紹介します．データ処理といっても，条件にあう列や行を抽出したりするだけです． ここでは，具体例として，ユニオンツールのmyBeatの心拍間隔データを使います．フ…

信号処理では，時系列に含まれる特定の周波数成分を取り出すために，ローパス (low pass)フィルタとか，ハイパス (high pass)フィルタとか，バンドパス (band pass)フィルタとかを使うことがあります．学生の皆さんには，デジタルフィルタの基礎の理論から実…

前回，フォルダ内にあるたくさんのファイルについて，フォルダ内にある全ファイル名を取得するRのコマンドとしてlist.files(...)の使い方を紹介しました．今回は，ファイル名に含まれる文字を使って，特定のファイルを見つける方法を説明します．今回，メイ…

コンピュータを使えば，大量のデータを一括処理できます．つまり，たくさんのデータファイルがあっても，それらを自動で読み込んで，分析できちゃいます．今回は，Rで，たくさんのデータファイルを処理するときのヒントです． 事前知識 Rでフォルダ内にある…

ひさびさに，共同研究でコサイナー解析 (cosinor analysis)を頼まれたので，ポイントをメモっておきます．コサイナー解析というのは，周期的に変化すると思われる時系列に，コサインの波形を当てはめてる分析法です．概日リズムの特徴付けに使われることがあ…

今回は，Savitzky-Golayフィルタの畳込み表現の係数 (下図のピンク図形)と周波数特性の話です． Savitzky-Golay平滑化フィルタ．元の時系列（上段灰色）は，下段の破線にノイズを加えたもの．下段の灰色破線はノイズを加える前の時系列．赤実線（上下両方）…

前回，振動している信号のエンベロープ (包絡線)をヒルベルト変換で求める方法を紹介しました．今回は，この方法がうまくいくための注意事項をいくつか書いておきます． ポイントは， 前処理として，周波数の狭いバンドパスをかけておく あるいは 前処理とし…

振動する時系列の中には，振動の振幅の変化に役立つ情報がある場合があります．例えば，呼吸信号（換気量や気流）は振動します．そのとき，心拍数に影響を与えるのは，振動の周期ではなく，振幅です．ということで，今回はその「振幅」の情報を取り出すため…

今回は，2変量の自己回帰過程 について最小2乗法でパラメタを推定し，AICで次数を決定します．お気楽な"vars"パッケージのVARなどは使いません．計算過程や式を確認できるように，Rスクリプトを書きました．Rスクリプトは一番下に掲載してあります．AICにつ…

今回は，時系列に対して自己回帰過程をフィットし，最適な次数を決定します．自己回帰過程の次数の推定には，赤池情報量規準 (AIC: Akaike Information Criterion)を使います．AICの導出をここで説明することは，私には難しいのでやりません．でも，計算過程…

今回は，コヒーレンス (coherence，あるいは，magnitude-squared coherence)の計算のメモです．コヒーレンスは，2つの時系列間の相関を周波数領域で見る方法です．2つの信号の周波数成分に，線形な入出力関係があれば1に近くなります．なければ0です．とはい…

観測された時系列に最小２乗法で自己回帰過程を当てはめる計算のメモです．どんな形の過程を仮定しても，基本的な考え方を理解しておけば，係数を求めることができます．楽をしたければ，arとか，varsパッケージなどを使えば良いと思います．しかし，お気楽…

前回に続き，ARFIMA(0, d, 0)を，無限次の自己回帰過程 (AR過程: autoregressive process)で表現する話です．今回は数値実験です． 1. 基本事項：long AR過程表現 2. 数値計算のポイント 3. 数値計算結果 Rスクリプト 前回の記事は，これです． chaos-kiyono…

今回はギザギザ，凸凹したデータの平滑化法の一つであるSavitzky-Golayフィルタの紹介です． Savitzky-Golayフィルタの考え方 (2次でスケール5の例)：部分区間に多項式をフィットして中央の値を平滑化された値として採用する．実用上は，畳み込み和を使って…

グレンジャー因果 (Granger causality)のGewekeのアプローチを少しずつ理解していきます．前回の説明の続きで，今回は数値実験をします． chaos-kiyono.hatenablog.com 1. 問題設定 単純な例から考えはじめないと，私には理解できません．一般的な表現を見て…

Rで時系列のパワースペクトルを推定したいとき，もっともお気楽なのは，"spectrum"というコマンドを使うことだと思います．今回は，心拍変動 (RR間隔時系列)のスペクトル解析で"spectrum"コマンドを使った実例を紹介します．ちなみに，例として使用したRR間…

Rのコマンド"spectrum"のオプション"spans"の設定で，パワースペクトルの推定結果がどうなるかのデモです．"spans"が何かは，前回の記事を参考にしてください． chaos-kiyono.hatenablog.com 時系列の生成については，以下の記事を参考にしてください． chao…