如何计算正弦曲线两端的面积

问题描述

给定这个数据集:

y<-c(-13,16,35,40,28,36,43,33,40,33,22,-5,-27,-31,-29,-25,-26,-31,-26,-24,-25,-29,-23,4)
t<-1:24

我的目标是计算两个区域.第一个区域将仅集成来自零线上方的曲线第一部分的数据.第二个区域将整合零线下方曲线第二部分的数据.

My goal is to calculate two areas. The first area would integrate only data from the first part of the curve found above the Zero line. The second area would integrate data from the second part of the curve found below the zero line.

首先，我想对这个数据拟合一个正弦波.使用这个优秀的答案:

First I would like to fit a sine wave to this data. Using this excellent answer:

https://stats.stackexchange.com/questions/60994/fit-a-sinusoidal-term-to-data

我能够拟合正弦波(我将使用具有二次谐波的周期，看起来更适合)

I was able to fit a sine wave (I will be using the periodic with second harmonic which looks to have a better fit)

ssp <- spectrum(y)
per <- 1/ssp$freq[ssp$spec==max(ssp$spec)]
reslm <- lm(y ~ sin(2*pi/per*t)+cos(2*pi/per*t))
summary(reslm)

rg <- diff(range(y))
plot(y~t,ylim=c(min(y)-0.1*rg,max(y)+0.1*rg))
lines(fitted(reslm)~t,col=4,lty=2)   # dashed blue line is sin fit

# including 2nd harmonic really improves the fit
reslm2 <- lm(y ~ sin(2*pi/per*t)+cos(2*pi/per*t)+sin(4*pi/per*t)+cos(4*pi/per*t))
summary(reslm2)
lines(fitted(reslm2)~t,col=3)    # solid green line is periodic with second harmonic
abline(h=0,lty=2)

接下来我想计算仅为正的曲线下面积，以及完全为负的曲线下面积.我很幸运地使用 Bolstad2 和 Mess 包中的 AUC 函数查看了类似的答案.但是我的数据点并没有整齐地落在零线上，我不知道如何将正弦函数分解为仅在零线上方和仅在零线下方的区域.

Next I would like to calculate the area under the curve that is only positive, as well as the area under the curve that is exclusively negative. I've had luck looking at similar answers using the AUC functions in the Bolstad2 and Mess packages. But my data points do not fall neatly on zero line, and I do not know how to break up the sine function into areas only above the Zero line and only below the Zero line.

推荐答案

第一件事.要进行精确计算，您需要使用 2 次谐波傅立叶的精确函数.其次，谐波函数的美妙之处在于它们是重复的.所以如果你想找到你的函数在哪里达到 0，你只需要将你的区间扩展到，这样你就可以确保找到 2 个以上的根.

First things first. To get an exact calculation, you will need to work with the exact function of the 2nd harmonic fourier. Secondly, the beauty of harmonics functions is that they are repetitive. So if you want to find where your function reaches 0, you merely need to expand your interval to so you can be sure to find more than 2 roots.

首先我们从回归模型中得到准确的函数

First we get the exact function from the regression model

fourierfnct <- function(t){
  fnct <- reslm2$coeff[1]+
    reslm2$coeff[2]*sin(2*pi/per*t)+
    reslm2$coeff[3]*cos(2*pi/per*t)+
    reslm2$coeff[4]*sin(4*pi/per*t)+
    reslm2$coeff[5]*cos(4*pi/per*t)
  return(fnct)
}

其次，您可以编写一个可以找到根的函数(函数为0).R 提供了一个 uniroot 函数，您可以使用它来查找循环中的多个根.

secondly,you can write a function which can find the roots (where the function is 0). R provides a uniroot function which you can use to find multiple roots in a loop.

manyroots <- function(f,inter,period){
  roots <- array(NA, inter)
  for(i in 1:(length(inter)-1)){
    roots[i] <- tryCatch({
      return_value <- uniroot(f,c(inter[i],inter[i+1]))$root
    }, error = function(err) {
      return_value <- -1
    })
  }
  retroots <- roots[-which(roots==-1)]
  return(retroots)
}

然后您只需计算根，并使用它们跨这些边界对函数进行积分.

then you simply calculate the roots, and use them to integrate the function across those boundaries.

roots <- manyroots(fourierfnct,seq(0,25),per)
integrate(fourierfnct, roots[1],roots[2])
#300.6378 with absolute error < 3.3e-12
integrate(fourierfnct, roots[2],roots[3])
#-284.6378 with absolute error < 3.2e-12

这篇关于如何计算正弦曲线两端的面积的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！