Datenanpassung an trigonometrische FunktionenHallo,
Ich habe Schwierigkeiten, gemessene Daten an eine trigonometrische Ksin(ax+b)+c-Funktion anzupassen. Das liegt in der Natur meines Prozesses und ich suche nach einer automatischen Möglichkeit, ihn per JSL-Code anzupassen.
Ich habe einen Code geschrieben, der die Daten normalisiert, sodass die Schwankung zwischen -1 und 1 liegt.
Am Ende möchte ich die Konstanten a und b erhalten (auch K und c, aber es ist viel einfacher)
beigefügtes Beispiel für meinen Datentyp.
BR,
Yoav
Automatisierung und SkripterstellungDatenexploration und -visualisierungQualitäts- und VerfahrenstechnikBetreff: Datenanpassung an trigonometrische FunktionenDu könntest einen gebrauchen Nulldurchgangsdetektor , was ein kleines Stück JSL wäre, um Ihre Daten zu durchsuchen und jedes Mal zu bemerken, wenn die Daten von negativ zu positiv wechseln. Das ist der Beginn eines Sinuswellenzyklus. Der Abstand zum nächsten Nulldurchgang ist die Periode.
Möglicherweise haben Sie Probleme mit verrauschten Daten, die sich scheinbar mehrmals überschneiden.
Der Nulldurchgang ist leichter zu erkennen als Spitzen.
Ich habe diese Idee verwendet Panne .
Betreff: Datenanpassung an trigonometrische FunktionenHallo, vielen Dank an alle für die Zeit, die ihr dafür aufgewendet habt. Das scheint ein toller Ansatz für meine Anwendung zu sein. Aber wie haben Sie die Periodenlänge „automatisch“ ermittelt? Das habe ich nicht verstanden..
Yoav
Betreff: Datenanpassung an trigonometrische FunktionenSuper, das ist genau die Intuition, die mir gefehlt hat.Stattdessen habe ich die Parameterwerte in der nichtlinearen Plattform manipuliert, bis das Diagramm den Daten nahe kam, aber zu wissen, was die Parameter darstellen (und aus den Daten gute erste Vermutungen zu ziehen), ist ein viel besserer Weg.
Betreff: Datenanpassung an trigonometrische FunktionenIch denke, es hilft, das zu verwenden, was ich als Standardform einer Sinusfunktion bezeichne: Y = K*sin(A*(x -
) + C.
In Ihren Daten:
K = 1/2 des Wertebereichs
A = 2*pi / Periode
B = Verschiebung nach rechts von 0 zum Startpunkt der Sinuswelle
C = Mittelwert Ihrer Y-Werte
Anfangsschätzungen für die nichtlineare Regression können mithilfe dieser Definitionen aus der Distributionsplattform und dem Graph Builder ermittelt werden:
K = 1/2 * (Max - Min) = 1/2*(47 - 17,6) = 14,7
A = 2 * 3,14 / 880 = 0,007
B = 662
C = 34,1


Die Nichtlinea r Rückschritt Bei diesem Ansatz konvergiert das Ion sehr schnell.Das Modell konvergierte in 5 Iterationen unter Verwendung dieser Werte. Siehe unten

Betreff: Datenanpassung an trigonometrische FunktionenHier erfahren Sie, wie Sie eine Schätzung erhalten A Parameter mithilfe der FFT-Funktion. Bei den Kommentaren zu Sekunden und Stunden wird davon ausgegangen, dass die Daten etwa eine Stunde lang etwa alle 4 Sekunden erfasst wurden. Die FFT kümmert sich nicht um die Einheiten. Die FFT klassiert die Frequenzen in den Daten so, dass das erste Bin nach dem DC-Offset-Bin die niedrigste Frequenz darstellt, die einen vollständigen Zyklus in den Eingabedaten erzeugen kann, das nächste Bin zwei vollständige Zyklen umfasst, ... das mittlere Bin ndata/ ist. 2 Zyklen. Die Behälter danach sind ein Spiegelbild des ersten.
dt = Open( "z:\example data.jmp" );
// by hand, I made a formula column to look at x-lag(x) and
// decided the first 79 and last 8 were not representative...
rowstart = 80;
rowstop = 937;
timeStart = dt:x[rowstart];
timeStop = dt:x[rowstop];
ydata = dt[rowstart :: rowstop, {y}]; // throw away a few with irregular time
{r, i} = FFT( {ydata} );
nBins = N Items( r );
sampleRate = (timeStop - timeStart) / (nBins - 1); // imagine there were only 2 bins...thus -1.
// the samplerate is ~4 seconds/sample...
// the highest frequency this can resolve is 1 cycle/ 8 seconds. (.125 Hz)
// frequencies higher than that will cause aliasing unless they are filtered out before this.
// the bin [nBins/2] represents the maximum frequency that can be resolved
lastBinFrequency = 1 / (2 * sampleRate); // ~ 0.125 Hz
// the first bin [1] is the DC, or offset from 0, component. We'll ignore it.
// the next bin [2] represents the frequency that has one cycle over the length of the data.
// the lowest non-DC frequency is one cycle per hour (the length of time the data was collected).
binWidth = 1 / (timeStop - timeStart); // ~0.0003Hz, or ~1 cycle per hour
// note: lastBinFrequency == (nbins-1)*binWidth/2 -- each bin is this wide, in frequency
// the remaining bins are a mirror of the first and ignored.
// combine real and imaginary bin components to get an amplitude
p = Sqrt( r ^ 2 + i ^ 2 );
dtGraph2 = As Table( ((rowstart :: rowstop) - rowstart)` || p );
dtGraph2 << Graph Builder(
Size( 1300, 500 ),
Show Control Panel( 0 ),
Show Legend( 0 ),
Variables( X( :Col1 ), Y( :Col2 ) ),
Elements( Points( X, Y, Legend( 5 ) ) ),
SendToReport(
Dispatch( {}, "Col1", ScaleBox, {Format( "Best", 12 ), Min( 0 ), Max( nBins / 2 ), Inc( 100 ), Minor Ticks( 4 )} ),
Dispatch( {}, "Col2", ScaleBox, {Format( "Best", 12 ), Min( 0 ), Max( 30000 ), Inc( 5000 ), Minor Ticks( 1 )} ),
Dispatch( {}, "graph title", TextEditBox, {Set Text( "FFT Bins" )} ),
Dispatch( {}, "X title", TextEditBox, {Set Text( "bin number" )} ),
Dispatch( {}, "Y title", TextEditBox, {Set Text( "size" )} )
)
);
// there is a stong DC component in this data (all values > 0) so bin 1 is the max.
If( Loc Max( p ) != 1,
Throw( "something odd" )
);
p[1] = 0; // clear the DC value
// next max after DC:
maxBinIndex = Loc Max( p ); // 5. the data repeats about (5-1)==4 times in this frame.
Write(
"\!nBiggest frequency component at ",
binWidth * (maxBinIndex - 1),
" Hz or ",
Round( 1 / (binWidth * (maxBinIndex - 1)) ),
" seconds per cycle"
);
// where does a==.007 come from?
secondsPerCycle = 1 / (binWidth * (maxBinIndex - 1));
totalSeconds = timeStop-timeStart;
cyclesInSample = totalSeconds/secondsPerCycle;//about 4
timerange = timeStop-timeStart;
a_from_fft = cyclesInSample*2*pi()/timerange; // .0073...
show(a_from_fft);
One strong response in 4th bin after DC binOne strong response in 4th bin after DC binEine starke Reaktion im 4. Bin nach DC-Bin
Größte Frequenzkomponente bei 0,00116448326055313 Hz oder 859 Sekunden pro Zyklus
a_from_fft = 0,007316664113164;
Die DC-Komponente hängt mit dem Mittelwert der Daten zusammen C Wert.
Betreff: Datenanpassung an trigonometrische FunktionenIch weiß nicht, ob Sie das wollen, aber Sie können das Problem mithilfe der nichtlinearen Plattform lösen.Im Anhang finden Sie Ihre Daten sowie eine Prädiktorspalte (unter Verwendung der zu schätzenden Parameter) und ein Skript, das diese Parameter auflöst.Ich gebe zu, dass es eine ganze Weile gedauert hat, bis ich erste Schätzungen über die Parameter erhalten habe, die die Konvergenz der Lösung ermöglichten (wahrscheinlich aufgrund meiner mangelnden Vertrautheit mit der Sinusfunktion – es ist lange her, dass ich mit solchen Funktionen gearbeitet habe).Aber konzeptionell wird die nichtlineare Plattform das erreichen, was Sie wollen.