Bien que j'ai lu ce post, je n'ai toujours aucune idée de comment l'appliquer à mes propres données et j'espère que quelqu'un pourra m'aider.
J'ai les données suivantes:
y <- c(11.622967, 12.006081, 11.760928, 12.246830, 12.052126, 12.346154, 12.039262, 12.362163, 12.009269, 11.260743, 10.950483, 10.522091, 9.346292, 7.014578, 6.981853, 7.197708, 7.035624, 6.785289, 7.134426, 8.338514, 8.723832, 10.276473, 10.602792, 11.031908, 11.364901, 11.687638, 11.947783, 12.228909, 11.918379, 12.343574, 12.046851, 12.316508, 12.147746, 12.136446, 11.744371, 8.317413, 8.790837, 10.139807, 7.019035, 7.541484, 7.199672, 9.090377, 7.532161, 8.156842, 9.329572, 9.991522, 10.036448, 10.797905)
t <- 18:65
Et maintenant je veux simplement rentrer dans une onde sinusoïdale
avec les quatre inconnues , , et à elle.
Le reste de mon code est le suivant
res <- nls(y ~ A*sin(omega*t+phi)+C, data=data.frame(t,y), start=list(A=1,omega=1,phi=1,C=1))
co <- coef(res)
fit <- function(x, a, b, c, d) {a*sin(b*x+c)+d}
# Plot result
plot(x=t, y=y)
curve(fit(x, a=co["A"], b=co["omega"], c=co["phi"], d=co["C"]), add=TRUE ,lwd=2, col="steelblue")
Mais le résultat est vraiment médiocre.
J'apprécierais beaucoup toute aide.
À votre santé.
r
regression
fitting
Pascal
la source
la source
Réponses:
Si vous voulez juste une bonne estimation de et que vous ne vous souciez pas beaucoup de son erreur standard:ω
(Un meilleur ajustement expliquerait peut-être encore les valeurs aberrantes de cette série, en réduisant leur influence.)
---
Si vous voulez avoir une idée de l'incertitude dans , vous pouvez utiliser la vraisemblance du profil ( pdf1 , pdf2 - les références pour obtenir des CI ou SE approximatifs à partir de la vraisemblance du profil ou de ses variantes ne sont pas difficiles à localiser)ω
(Alternativement, vous pouvez alimenter ces estimations dans nls ... et les démarrer déjà convergées.)
la source
lm(y~sin(2*pi*t)+cos(2*pi*t)
mais cela n'a pas fonctionné (lecos
terme était toujours 1). Juste par curiosité: que font les deux premières raies (je sais quispectrum
estime la densité spectrale)?2*pi*t
spec
dans TSA peut être meilleure (elle semble avoir plus d'options, dont l'une peut parfois être importante), mais dans ce cas, le pic principal était exactement au même endroit que pour,spectrum
donc je n'ai pas pris la peine.reslm
enreslm <- lm(y ~ cos(2*pi/per*t)+tan(2*pi/per*t))
mais cela ne semble pas correct. des indices?Quand j'ai mis cela dans
nls
lastart
liste de, j'ai obtenu une courbe qui était beaucoup plus raisonnable, bien qu'elle ait encore des biais systématiques.En fonction de votre objectif avec cet ensemble de données, vous pouvez essayer d'améliorer l'ajustement en ajoutant des termes supplémentaires ou en utilisant une approche non paramétrique comme un processus gaussien avec un noyau périodique.
Choisir automatiquement une valeur de départ
Si vous souhaitez choisir la fréquence dominante, vous pouvez utiliser une transformée de Fourier rapide (FFT). C'est loin de mon domaine d'expertise, donc je vais laisser d'autres personnes remplir les détails s'ils le souhaitent (en particulier sur les étapes 2 et 3), mais le
R
code ci-dessous devrait fonctionner.Vous pouvez également tracer
abs(truncated.fft)
pour voir s'il existe d'autres fréquences importantes, mais vous devrez jouer un peu avec la mise à l'échelle de l'axe des x.De plus, je crois que @Glen_b a raison que le problème est convexe une fois que vous connaissez les oméga (ou peut-être avez-vous aussi besoin de connaître le phi? Je ne suis pas sûr). Dans tous les cas, connaître les valeurs de départ pour les autres paramètres ne devrait pas être aussi important que pour les oméga s'ils sont dans le bon stade. Vous pourriez probablement obtenir des estimations décentes des autres paramètres de la FFT, mais je ne sais pas comment cela fonctionnerait.
la source
foo.bar
. Cela est dû à la façon dont R spécifie les méthodes pour les classes .Comme alternative à ce qui a déjà été dit, il peut être intéressant de noter qu'un modèle AR (2) de la classe des modèles ARIMA peut être utilisé pour générer des prévisions avec un modèle d'onde sinusoïdale.
Panratz (1991) nous apprend ce qui suit sur les cycles stochastiques:
Pour voir si un tel modèle pouvait être ajusté aux données, j'ai utilisé la
auto.arima()
fonction du package de prévision pour savoir s'il suggérerait un modèle AR (2). Il s'avère que laauto.arima()
fonction suggère un modèle ARMA (2,2); pas un pur modèle AR (2), mais c'est OK. C'est OK car un modèle ARMA (2,2) contient un composant AR (2), donc la même règle (sur les cycles stochastiques) s'applique. Autrement dit, nous pouvons toujours vérifier la condition susmentionnée pour voir si des prévisions d'ondes sinusoïdales seront produites.Les résultats de
auto.arima(y)
sont présentés ci-dessous.Le graphique ci-dessous montre la série d'origine, y, l'ajustement du modèle ARMA (2,2) et 14 prévisions hors échantillon. Comme on peut le voir, les prévisions hors échantillon suivent un schéma sinusoïdal.
Gardez à l'esprit deux choses. 1) Il s'agit simplement d'une analyse très rapide (à l'aide d'un outil automatisé) et un traitement approprié impliquerait de suivre la méthodologie de Box-Jenkins. 2) Les prévisions ARIMA sont bonnes pour les prévisions à court terme, vous pouvez donc trouver que les prévisions à long terme des modèles dans les réponses de @David J. Harris et @Glen_b sont plus fiables.
Enfin, j'espère que c'est un bon ajout à certaines réponses déjà très instructives.
Référence : Prévision avec des modèles de régression dynamique: Alan Pankratz, 1991, (John Wiley and Sons, New York), ISBN 0-471-61528-5
la source
Les méthodes actuelles pour adapter une courbe sin à un ensemble de données donné nécessitent une première estimation des paramètres, suivie d'un processus interactif. Il s'agit d'un problème de régression non linéaire. Une autre méthode consiste à transformer la régression non linéaire en une régression linéaire grâce à une équation intégrale pratique. Ensuite, il n'y a pas besoin de conjecture initiale et pas besoin de processus itératif: l'ajustement est directement obtenu. Dans le cas de la fonction y = a + r * sin (w * x + phi) ou y = a + b * sin (w * x) + c * cos (w * x), voir pages 35-36 du document "Régression sinusoidale" publiée sur Scribd: http://www.scribd.com/JJacquelin/documents Dans le cas de la fonction y = a + p * x + r * sin (w * x + phi): pages 49-51 du chapitre "Régressions mixtes linéaires et sinusoïdales". Dans le cas de fonctions plus compliquées, le processus général est expliqué dans le chapitre "Régression sinusoïdale généralisée" pages 54-61, suivi d'un exemple numérique y = r * sin (w * x + phi) + (b / x) + c * ln (x), pages 62-63
la source
Si vous connaissez le point le plus bas et le plus haut de vos données d'aspect cosinus, vous pouvez utiliser cette fonction simple pour calculer tous les coefficients cosinus:
Ci-dessous, il est utilisé pour simuler la variation de la température tout au long de la journée avec une fonction cosinus, en entrant les heures et les valeurs de température pour l'heure la plus basse et la plus chaude:
La sortie est ci-dessous:
la source
Une autre option utilise la fonction générique optim ou nls. J'ai essayé les deux aucun d'eux n'est complètement robuste
Les fonctions suivantes prennent les données en y et calculent les paramètres.
l'utilisation est la suivante:
Le code suivant compare les données
la source