如果您真的不知道周期性是什么，那么最好的方法可能是找到对应于频谱密度最大值的频率。但是，低频的频谱会受到趋势的影响，因此您需要先对序列进行去趋势。以下 R 函数应该为大多数系列完成工作。它远非完美，但我已经在几十个示例上对其进行了测试，它似乎工作正常。对于没有强周期性的数据，它将返回 1，否则返回周期长度。

更新：功能的第 2 版。这要快得多，而且似乎更健壮。

    find.freq <- function(x)
    {
        n <- length(x)
        spec <- spec.ar(c(x),plot=FALSE)
        if(max(spec$spec)>10) # Arbitrary threshold chosen by trial and error.
        {
            period <- round(1/spec$freq[which.max(spec$spec)])
            if(period==Inf) # Find next local maximum
            {
                j <- which(diff(spec$spec)>0)
                if(length(j)>0)
                {
                    nextmax <- j[1] + which.max(spec$spec[j[1]:500])
                    period <- round(1/spec$freq[nextmax])
                }
                else
                    period <- 1
            }
        }
        else
            period <- 1
        return(period)
    }

如果您期望该过程是静止的——周期性/季节性不会随着时间而改变——那么像卡方周期图（参见 Sokolove 和 Bushell，1978 年）之类的东西可能是一个不错的选择。它通常用于分析昼夜节律数据，其中可能包含大量噪声，但预计具有非常稳定的周期性。

这种方法不对波形的形状做任何假设（除了它在各个周期之间是一致的），但确实要求任何噪声都具有恒定的平均值并且与信号不相关。

chisq.pd <- function(x, min.period, max.period, alpha) {
N <- length(x)
variances = NULL
periods = seq(min.period, max.period)
rowlist = NULL
for(lc in periods){
    ncol = lc
    nrow = floor(N/ncol)
    rowlist = c(rowlist, nrow)
    x.trunc = x[1:(ncol*nrow)]
    x.reshape = t(array(x.trunc, c(ncol, nrow)))
    variances = c(variances, var(colMeans(x.reshape)))
}
Qp = (rowlist * periods * variances) / var(x)
df = periods - 1
pvals = 1-pchisq(Qp, df)
pass.periods = periods[pvals<alpha]
pass.pvals = pvals[pvals<alpha]
#return(cbind(pass.periods, pass.pvals))
return(cbind(periods[pvals==min(pvals)], pvals[pvals==min(pvals)]))
}

x = cos( (2*pi/37) * (1:1000))+rnorm(1000)
chisq.pd(x, 2, 72, .05)

最后两行只是一个例子，表明它可以识别纯三角函数的周期，即使有很多加性噪声。

如所写，alpha调用中的最后一个参数 ( ) 是多余的，该函数只是返回它可以找到的“最佳”周期；取消注释第一个return语句并注释掉第二个语句，使其返回该级别所有重要期间的列表alpha。

此功能不会进行任何类型的完整性检查以确保您已放入可识别的句点，它（可以）也不适用于小数句点，如果您决定，也没有内置任何类型的多重比较控制查看多个时期。但除此之外，它应该是相当稳健的。

你可能想更清楚地定义你想要的东西（如果不是在这里，对你自己）。如果您要寻找的是噪声数据中包含的最具统计意义的静止期，那么基本上有两条路线可供选择：

1) 计算稳健的自相关估计，取最大系数
2) 计算稳健的功率谱密度估计，取频谱的最大值

#2 的问题在于，对于任何嘈杂的时间序列，您将在低频下获得大量功率，从而难以区分。有一些技术可以解决这个问题（即预白化，然后估计 PSD），但如果您的数据的真实周期足够长，那么自动检测将是不确定的。

您最好的选择可能是实现一个健壮的自相关例程，例如Maronna、Martin 和 Yohai的Robust Statistics - Theory and Methods的第 8.6、8.7 章。在 Google 上搜索“robust durbin-levinson”也会产生一些结果。

如果您只是在寻找一个简单的答案，我不确定是否存在。时间序列中的周期检测可能很复杂，并且要求可以执行魔术的自动化例程可能太多了。

您可以使用 DSP 理论中的希尔伯特变换来测量数据的瞬时频率。网站http://ta-lib.org/有开源代码，用于衡量金融数据的主导周期；相关函数称为HT_DCPERIOD；您也许可以使用它或根据您的目的调整代码。