0
$\begingroup$

I have the following picture: enter image description here

I want to use R language for estimating the potential wave. I am looking for answers for the following questions:

  1. Is this the wave ?
  2. The goodness of fit my wave describes the data? (Like $R^2$)
  3. Does the wave rise, fall or stay constant ?
  4. What are the local and absolute extreme points of the wave?

  • Which model is best to make such a measurement ? Is a non-linear regression a good option
  • Which library in R do you recommend?
  • Which indicators are useful measures for my problem (e.g. p-values)?
$\endgroup$
5
$\begingroup$

That is a very nice question to deal with the topic of regression.

First and foremost it is nescessary to know what you want to know about your curve and imho you are asking if you could extract some basic information about the sinusoid (a mathematical curve that describes a smooth repetitive oscillation). The main interest of yours seems to be wether or not the slope of this sinusoid is positive, negative, or 0.

Lets produce some data first:

library(tidyverse)

##
# the linear case
##
someData <- data.frame(x = runif(200, -4*pi, 4*pi)) %>%
  mutate(y = 5 - 3*x + 10*sin(x) + rnorm(n(), sd = 3))

ggplot(someData, aes(x,y)) +
  geom_point()

wavy

Seems wavyy enough to me ;-) - it has a linear component, an intercept, and an oscilating behaviour along the linear component.

Here are some approaches to think about:


Linear regression

You are dealing with linear regression if you assume your coefficients $\vec{\beta}$ to be linear (and it does not matter that the variables themselfes are not), i.e. you are trying to find a solution to

$$ y = \beta_0 + \sum\limits_i \beta_i*f(x_i) $$

Lets consider (what we already know) that our $y$ is the function of

$$ y = \beta_0 + \beta_1 x + \beta_2 sin(x) $$

Here our $\vec{\beta}$s are linear and thus we can fit a linear model:

myLinearModel <- lm(y ~ x + sin(x), someData)
# Coefficients:
#(Intercept)            x       sin(x)  
#      4.795       -2.944        9.864  
summary(myLinearModel)
#Coefficients:
#            Estimate Std. Error t value Pr(>|t|)    
#(Intercept)  4.79516    0.20465   23.43   <2e-16 ***
#x           -2.94381    0.03051  -96.49   <2e-16 ***
#sin(x)       9.86405    0.30394   32.45   <2e-16 ***

A very good fit, indeed (as related to question 2). It already answers your question 3 about what direction the wave is going (-2.994), i.e. it has a negative slope in the linear component. But for now you cannot say anything about contraction or expansion.

We can also plot the fit with the data:

ggplot(someData, aes(x,y)) +
  geom_point() + 
  geom_smooth(method = "lm", formula = y ~ x + sin(x))

enter image description here


Non-linear regression

Lets construct some non-linear data

someData <- data.frame(x = runif(500, -4*pi, 4*pi)) %>%
  mutate(y = 5 - 3*x + 5*sin(x) + 5 * cos(2*x + 2) + rnorm(n(), sd = 1))

ggplot(someData, aes(x,y)) +
  geom_point()

enter image description here

We still can try and fit a linear model

myLinearModel <- lm(y ~ x + sin(x) + cos(x), someData)
#Coefficients:
#Coefficients:
#(Intercept)            x       sin(x)       cos(x)  
#    5.07889     -2.89829      5.03452     -0.08321
summary(myLinearModel)
#Coefficients:
#            Estimate Std. Error  t value Pr(>|t|)    
#(Intercept)  5.07889    0.16187   31.376   <2e-16 ***
#x           -2.89829    0.02152 -134.687   <2e-16 ***
#sin(x)       5.03452    0.22566   22.311   <2e-16 ***
#cos(x)      -0.08321    0.23439   -0.355    0.723   

Still, a good fit, but its seems the cos component is negligible.

Now, you can still have your question about the slope of the linear component answered. Yay!

But if you would like to inter- or extrapolate some new points you can build different models. E.g. loess

ggplot(someData, aes(x,y)) +
  geom_point() + 
  geom_smooth(method = "lm", formula = y ~ x + sin(x) + cos(x)) +
  geom_smooth(method = "loess", col = "red", span = 0.1) 

enter image description here

With a small span you actually can very well use it for interpolation - but the low the span the greater the overfitting!

Another form of non-linear models can be used in general optimization procedures, e.g.:

 nonLinear <- nls(y ~ i + x + a * sin(x) + b * cos(c * x), data = someData,
             start = structure(as.list(runif(7, -10, 10)), names = c("a", "b", "c", "i")))
nonLinear
summary(nonLinear)

##
# a complex formula that probably will not find an optimum
##
nonLinear <- nls(y ~ i + x + a * sin(b * x + e1) + c * cos(x*d + e2), data = someData,
                 start = structure(as.list(runif(7, -10, 10)), names = c("a", "b", "c", "d", "e1", "e2", "i")))
summary(nonLinear)

I will not go into detail here; there is soooo much to say and think about it and it all depends on your data. And I will not dare to start with Decision Trees or Neural Networks for regression. (Because if you have a hammer everything seems like a nail)

Fourier Transform

What else? Well, how about some stuff from signal processing where wavyy functions are everywhere!

Using fft you could segment your function into a series of sinusoids like here:

fftin <- arrange(someData, x)$y
fftres <- fft(fftin)
fftres
barplot(Mod(fftres))

enter image description here

In this frequency plot we can see many small frequencies. Lets remove them and see how our new function fits the data.

fftcut <- fftres
fftcut[(Mod(fftres) < 100)] <- 0

plot(Re(fft(fft(fftin), inverse = TRUE)/length(fftres)), type="p")
lines(Re(fft(fftcut, inverse = TRUE)/length(fftres)), type="l", col="red")

enter image description here

We can even cut out more and still have a smooth curve.

fftcut <- fftres
fftcut[(Mod(fftres) < 1000)] <- 0

plot(Re(fft(fft(fftin), inverse = TRUE)/length(fftres)), type="p")
lines(Re(fft(fftcut, inverse = TRUE)/length(fftres)), type="l", col="red")

enter image description here


Enough of magic :-). To answer your questions you can work with a linear model. Even expansion and contraction can be modelled, albeit not without addiotional data wrangling. Peaks of your functions can be derived on a piece of paper after you compute your coefficients in a linear model.

$\endgroup$
  • $\begingroup$ can you explain more about the fft? where what is the model (with Fourier basis) used in fft. why fft and inverse fft will work? $\endgroup$ – Haitao Du Oct 19 '16 at 14:36
  • $\begingroup$ I don't have strong background knowledge about FFT but this resource di.fc.ul.pt/~jpn/r/fourier/fourier.html + di.fc.ul.pt/~jpn/r/fourier/fourier2.html has much practical info. Basically all this guy's markdowns are awesome ;-) $\endgroup$ – Drey Oct 19 '16 at 15:11
  • $\begingroup$ Hi, sorry I have no deep knowledge in R or Mathematic, can you explain please: 1) how you will calculate extreme points ?, I think if I find extreme points will be enough to make decision about type of wave 2) Combination like this "x + sin(x) + cos+ runif" is very important for fitting model, but I have no idea how it working :( 3) the formula "lm(y ~ x + sin(x), someData)" is linear regression ??, this is curve, or.. $\endgroup$ – Max Usanin Oct 19 '16 at 16:03
  • $\begingroup$ and one more 4) "overfitting" - how can it be in this case ? $\endgroup$ – Max Usanin Oct 19 '16 at 16:10
  • $\begingroup$ @MaxUsanin 1) you will need to find critical points of trigonometric function, i.e. $f'(x) = \beta_1 - \beta_2 cos(x) = 0$ 2) Simply imagine the first formula for more $i$ and replace the function accordingly. rnorm is simply the noise in the data. 3) yes this is linear regression and the solution (i.e. the coefficients) describe the curve. 4) Overfitting relates to inter- and extrapolation. Try to add more noise and you will see that loess will start modeling noise as well. $\endgroup$ – Drey Oct 20 '16 at 7:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.