問題タブ [periodicity]

MATLABで微分方程式を解きたい

ここでdt、MyDiff、tSpanおよびX0は、コードの前半で定義されています。問題は、離散化ステップdtが非常に小さいのに、合計シミュレーション時間tSpan(end)が非常に長いことです。さらに、解はある時刻 T から周期 で周期的になると予想されPます。ここで、期間Pはアプリオリにわかっていますが、時間Tはわかりません。

私がやりたいことは、X計算時間を節約するために、解が周期的になったときに ode113-solver を自動的に停止することです。これをどのように行うべきかについての考えをいただければ幸いです。

今までの私の考え：

1. 最初の副問題は、周期性が発見されたときに Matlab ソルバーを停止する方法です。MATLAB には、イベント関数によって ode113-solver を停止するオプションが含まれています。

   odeopts = odeset('MaxStep',dt,'Events',MyEventFcn);

ただし、MyEventFcn関数でtありX、現在の時間ステップである必要があります。この情報から周期性を判断することは不可能のようです。X唯一の方法は、予想される期間の 2 倍までの以前の値を含むグローバル パラメータを使用することPです。しかし、これは私には非常に洗練されておらず、非効率的であるように思われ、より良い方法があることを望んでいました.

2. X2 番目の問題は、解が period で周期的であることをどのように判断するかPです。現時点では、自己相関関数にしきい値を使用するのが最善の方法だと考えていますがxcorr(X)、詳細についてはまだ確信が持てないため、ここでの提案も役に立ちます。

編集: 以前の考えにほぼ従って、グローバル変数を使用して X 値を追跡し、xcorr を使用して周期性を検出するソリューションを実装しました。しばらく考えてみると、実際にはこれが最も簡単な方法のように思えます。

以下は、多くの列を含む csv ファイルのスナップショットです (実際のファイルは非常に大きくなります)。各列は、「システム状態」のタイムスタンプ値に対応しています。私が見つけようとしているのは、これらの列の発生の期間です。つまり、システム状態がそれ自体を繰り返し、どの期間で発生するかです。私は python で fft2 を探していましたが、今ではピリオド抽出を理解していません。私はフーリエ変換に慣れておらず、それについての以前の知識もないので、助けてください。

1 つの行列が 1 つの列で表されます。最初の 2 列は、マトリックス セルの識別用です。ほとんどの値はゼロですが、すべてではありません。

私のプログラムのアルゴリズムステップ

周期成分分析を使用して、EEG チャネルの周期成分を検出したいと考えています。アルゴリズムはあまり知られていないので、手順を簡単に説明します。

周期性誤差関数:を最小化します。error(τ)= E[s(t+τ) -s(t)] / E[s(t)^2]ここτで、 はタイムラグです。行列A(τ), , を定義するとC、error(τ)=(w'*A(τ)*w)/(w'*C*w) = Rayleigh(A(τ), C, w). したがって、最小値は、特定のタイム ラグ τ に対する( Α(τ), ) の最小の一般化された固有値です。C

私はタイムラグの範囲で実行piCAし、通常、最小化された誤差関数の極小値で周期成分の周期を検出します (最小の一般化された固有値を使用)。 [E,W]=piCA( Xeeg, [minlag maxlag], 'pre', [1 1 1]); これは、変換行列Wと、Eすべてのペア ( A[τ], C) の一般化された固有値の行列を与えます。

最小の固有値で計算された関数ではなく、2 番目に小さい固有値で計算された関数の極小値で、正しい周期成分を検出することは可能ですか??

以下の画像を参照してください。

緑色の線は明らかにすべてのタイム ラグの最小誤差関数ですが、局所的な最小値はありません!!! したがって、結果は、2 番目に小さい固有値を使用して計算される青い線関数からのみ発生します。:/

半年ごと (半年ごと) のデータ ポイントを含む時系列があります。

「頻度 = 2」が実際の期間をはるかに超える非常に奇妙な時系列オブジェクトを返すため、ts() 関数はそれを処理できないようです。

Rでこの種の時系列オブジェクトの時系列分析を行う方法はありますか?

編集:ここに例があります:

私は期待していました：