Engee 文档
Notebook

趋势和谐波成分（傅立叶分解）

时间序列分析是研究各种过程的重要工具，无论是环境测量、设备的技术特性还是生物体的生命体征。

时间序列处理的关键阶段之一是识别隐藏在测量数据中的模式、趋势和周期成分。有了这些成分，就能更好地了解所研究现象的动态，并做出明智的预测。

本工作致力于通过光谱分析方法研究汽车尾气中所含一氧化碳（CO）浓度的时间相关性。为了完成这项任务，我们开发了一种特殊算法，其中包括几个连续步骤。让我们一一列举。

1.数据解析：从文件中读取原始数据，并将其转换为便于进一步处理的格式。

2.添加时间刻度：计算每个观测点相对于测量开始的第二个刻度。

3.趋势提取：通过对滑动窗口中的数值进行平均，然后用三次多项式进行近似，从而确定一氧化碳水平的长期趋势。

4.频谱分析：应用快速傅立叶变换（FFT）来识别重要的波动频率。

5.信号重建：去除趋势后，重建构成过程结构的主要谐波成分。

6.最终可视化：将重建后的信号与原始数据进行比较，评估所建议方法的质量。

在本文中，我们将以实际实验数据为例，展示分析方法识别时间序列特征的可能性。

算法实施说明

算法的主要阶段

1.数据准备

为了准备数据，我们使用了函数parse_line ，该函数可以将字符串分割成独立的数值，并将其转换为相应的数据类型。解析后，数据被加载到DataFrame 类型的表中。

然后添加秒格式的时间戳，时间戳是相对于第一次测量计算的。

2. 趋势检测

使用两种趋势检测方法：

  • 固定大小滑动窗口中的平均值、
  • 用三次多项式近似。

第一种方法简单方便，可以直观地估计总体趋势，而第二种方法可以精确地建立参数模型。

3.频谱分析

消除趋势后，应用快速傅立叶变换（FFT）。因此，可以确定特定时间序列的主要频率特征。在对获得的频谱进行分析时，可以识别出三个最大的分量（谐波）。

4.信号重构

根据选定的谐波，重建近似信号，并将其叠加到原始曲线上，然后根据趋势评估重建的准确性。

In [ ]: 
# Функция для парсинга одной строки
function parse_line(line)
    # Заменяем запятые на точки для корректного парсинга чисел
    line = replace(line, "," => ".")
    # Разделяем по пробелам/табам и фильтруем пустые элементы
    elements = filter(!isempty, split(line, r"\s+"))
    # Преобразуем элементы в нужные типы
    co = parse(Float64, elements[1])
    ch = parse(Float64, elements[2])
    co2 = parse(Float64, elements[3])
    o2 = parse(Float64, elements[4])
    lamb = parse(Float64, elements[5])
    n = parse(Int, elements[6])
    nox = parse(Int, elements[7])
    tmas = parse(Float64, elements[8])
    time = Time(elements[9])
    return (CO=co, CH=ch, CO2=co2, O2=o2, lamb=lamb, n=n, NOx=nox, Tmas=tmas, Time=time)
end
Out[0]:
parse_line (generic function with 1 method)
In [ ]: 
using DataFrames
# Чтение данных
data_lines = readlines("$(@__DIR__)/data.txt")[2:end]  # Пропускаем заголовок
parsed_data = [parse_line(line) for line in data_lines] # Парсим все строки
df = DataFrame(parsed_data) # Преобразуем в DataFrame
println(first(df, 5)) # Выводим первые 5 строк для проверки
println()
using Dates, FFTW, Statistics
# Добавим колонку с временем в секундах
start_time = df.Time[1]
df[!, :Seconds] = [Dates.value(t - start_time) / 1e9 for t in df.Time]  # переводим наносекунды в секунды
println(first(df, 5)) # Выводим первые 5 строк для проверки

5×9 DataFrame
 Row │ CO       CH       CO2      O2       lamb     n      NOx    Tmas     Time
     │ Float64  Float64  Float64  Float64  Float64  Int64  Int64  Float64  Time
─────┼──────────────────────────────────────────────────────────────────────────────
   1 │    0.64    113.0    13.83     0.75    1.012    820     13      0.0  00:00:00
   2 │    0.64    112.0    13.83     0.75    1.012    770    185      0.0  15:56:25
   3 │    0.64    112.0    13.83     0.75    1.012    790      0      0.0  15:56:26
   4 │    0.64    112.0    13.83     0.75    1.012    800      0      0.0  15:56:27
   5 │    0.64    112.0    13.84     0.74    1.012    990      0      0.0  15:56:28

5×10 DataFrame
 Row │ CO       CH       CO2      O2       lamb     n      NOx    Tmas     Time      Seconds
     │ Float64  Float64  Float64  Float64  Float64  Int64  Int64  Float64  Time      Float64
─────┼───────────────────────────────────────────────────────────────────────────────────────
   1 │    0.64    113.0    13.83     0.75    1.012    820     13      0.0  00:00:00      0.0
   2 │    0.64    112.0    13.83     0.75    1.012    770    185      0.0  15:56:25  57385.0
   3 │    0.64    112.0    13.83     0.75    1.012    790      0      0.0  15:56:26  57386.0
   4 │    0.64    112.0    13.83     0.75    1.012    800      0      0.0  15:56:27  57387.0
   5 │    0.64    112.0    13.84     0.74    1.012    990      0      0.0  15:56:28  57388.0

接下来，我们选择要分析的变量。让我们以 CO 为变量（可以用任何其他变量代替）：

In [ ]: 
i = 1 # Номер столбца для анализа
y = df[:,i]
column_names = names(df) # Получаем имена столбцов
println("Имя столбца для анализа: $(column_names[i])")
t = df.Seconds.-(15*60*60+56*60); t[1] = 0;
println("t: $(t[1:5])") # Выводим первые 5 строк для проверки

Имя столбца для анализа: CO
t: [0.0, 25.0, 26.0, 27.0, 28.0]

接下来，我们绘制趋势图。我们使用移动平均线来突出趋势：

In [ ]: 
window_size = 15  # размер окна для сглаживания
trend = [mean(y[max(1, i-window_size):min(end, i+window_size)]) for i in 1:length(y)]

# График исходных данных и тренда
plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO", legend=:topleft)
plot!(t, trend, label="Тренд", linewidth=2)
Out[0]:

为了获得更准确的趋势，我们可以使用多项式近似值来代替移动平均线：

In [ ]: 
using Polynomials
trend_coef = fit(t, y, 3)  # полином третьей степени
trend = trend_coef.(t)

# График исходных данных и тренда
plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO", legend=:topleft)
plot!(t, trend, label="Тренд", linewidth=2)
Out[0]:

让我们对谐波成分进行傅立叶分解。

In [ ]: 
# Удаляем тренд для анализа периодических составляющих
detrended = y .- trend

# Вычисляем БПФ
n = length(detrended)
freq = fftfreq(n, 1/mean(diff(t)))  # частоты в Гц
fft_vals = fft(detrended)

# Амплитуды и фазы
amplitude = abs.(fft_vals) ./ n * 2
phase = angle.(fft_vals)

# Оставляем только положительные частоты (симметрия БПФ)
pos_freq = freq[1:div(n,2)+1]
pos_ampl = amplitude[1:div(n,2)+1]

# График амплитудного спектра
plot(pos_freq, pos_ampl, xlabel="Частота (Гц)", ylabel="Амплитуда", 
     title="Амплитудный спектр", legend=false)
Out[0]:

重建主要谐波。让我们选择一些最重要的频率，恢复它们的贡献：

In [ ]: 
# Находим 3 наиболее значимые частоты
top_n = 3
sorted_idx = sortperm(pos_ampl, rev=true)
top_freq = pos_freq[sorted_idx[1:top_n]]
top_ampl = pos_ampl[sorted_idx[1:top_n]]
top_phase = phase[sorted_idx[1:top_n]]

# Восстанавливаем гармоники
harmonics = zeros(n)
for (f, a, ϕ) in zip(top_freq, top_ampl, top_phase)
    harmonics .+= a .* cos.(2π * f * t .+ ϕ)
end

# Итоговый график
plot(t, detrended, label="Детрендированные данные", xlabel="Время (с)", ylabel="CO")
plot!(t, harmonics, label="Основные гармоники", linewidth=2)
Out[0]:

全面分解：趋势和谐波。

In [ ]: 
plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO")
plot!(t, trend, label="Тренд", linewidth=2)
plot!(t, trend .+ harmonics, label="Тренд + гармоники", linewidth=2, linestyle=:dash)
Out[0]:

结论

在分析过程中，可以确定所研究时间序列行为的以下重要方面。

  • 确定了一氧化碳水平的长期波动，其明显趋势由车辆行驶和外部环境因素决定。
  • 确定了一氧化碳浓度周期性变化的结构，该结构通过三个主要振荡频率的分配来表示。
  • 对原始信号的重建质量进行了评估，结果表明去趋势数据与主要谐波之间具有良好的相关性。

因此，事实证明，所开发的方法对于分析物质浓度的时间依赖性等信号非常有效，研究人员甚至可以从真实实验的嘈杂数据中提取有用信息。

所提出的算法可以很容易地应用于分析其他研究领域的类似时间序列，在这些领域，将数据分离为趋势和高频波动非常重要。