SciPy周期图与AstroPy Lomb-Scargle周期图低频段结果差异

我分别使用SciPy的periodogram和AstroPy的Lomb-Scargle周期图计算数据的周期图，结果显示除低频段（接近最小频率）外，二者在其他频段均匹配，这是数值模拟的结果。

基于观测数据，我预期在0频率附近存在强信号，因此SciPy周期图的结果看起来更符合物理合理性。我尚未找到差异产生的原因及使二者结果一致的方法，恳请各位提供相关见解。

计算结果图像 #

SciPy周期图结果：

Lomb-Scargle周期图结果：

复现代码 #

from astropy.timeseries import LombScargle
import numpy as np
import pandas as pd
from scipy import signal
import requests 
import matplotlib.pyplot as plt

def plot_periodogram(x,y,N_freq,min_freq,max_freq,height_threshold,periodogram_type): 
    fig, ax = plt.subplots(figsize=(12,8))

    if periodogram_type == 'periodogram':
        dx = np.mean(np.diff(x))  # 假设x是均匀采样的
        fs = 1 / dx

        freq, power_periodogram = signal.periodogram(y,fs,scaling="spectrum",nfft=N_freq,
                                                     return_onesided=True,detrend='constant')
        power_max = power_periodogram[~np.isnan(power_periodogram)].max()
        
        plt.plot(freq, power_periodogram/power_max,linestyle="solid",color="black",linewidth=2)
        
        filename = "PowerSpectrum"
        
    else:
        
        freq = np.linspace(min_freq,max_freq,N_freq)
        ls= LombScargle(x, y,normalization='psd',nterms=1)
        power_periodogram= ls.power(freq)
                
        power_max = power_periodogram[~np.isnan(power_periodogram)].max()
        
        false_alarm_probabilities = [0.01,0.05]
        periodogram_peak_height= ls.false_alarm_level(false_alarm_probabilities,minimum_frequency=min_freq, 
                                                      maximum_frequency=max_freq,method='bootstrap')
        
        filename = "PowerSpectrum_LombScargle"
        plt.plot(freq, power_periodogram/power_max,linestyle="solid",color="black",linewidth=2)
        plt.axhline(y=periodogram_peak_height[0]/power_max, color='black', linestyle='--')
        plt.axhline(y=periodogram_peak_height[1]/power_max, color='black', linestyle='-')

    peaks_index, properties = signal.find_peaks(power_periodogram/power_max, height=height_threshold)    
    peak_values = properties['peak_heights']
    peak_power_freq = freq[peaks_index]

    for i in range(len(peak_power_freq)):
        plt.axvline(x = peak_power_freq[i],color = 'red',linestyle='--')
        ax.text(peak_power_freq[i]+0.05, 0.95, str(round(peak_power_freq[i],2)), color='red',ha='left', va='top', rotation=0,transform=ax.get_xaxis_transform())

    fig.patch.set_alpha(1)   
    plt.ylabel('Spectral Power',fontsize=20)
    plt.xlabel('Spatial Frequency', fontsize=20)
    plt.grid(True)
    plt.xlim(left=min_freq,right=max_freq)
   
    plt.xticks(fontsize=20)
    plt.yticks(fontsize=20)
    plt.savefig(filename,bbox_inches='tight')
    plt.show()

# CSV文件地址
url = 'https://pastebin.com/raw/uFi8WPvJ'

# 获取数据
response = requests.get(url)

if response.status_code == 200:
    data = response.text
    
    # 保存为CSV文件
    with open('data.csv', 'w') as f:
        f.write(data)
        
df =pd.read_csv('data.csv',sep=',',comment='%', names=['x', 'Bphi','r','theta'])
x = df['x'].values
y = df['Bphi'].values

# 删除含NaN的元素
indices = np.logical_not(np.logical_or(np.isnan(x), np.isnan(y)))
x = x[indices]
y = y[indices]

# 去均值
y = y - np.mean(y)

N_freq = 10000

min_freq = 0.001; 
max_freq = 4.0
height_threshold =0.7

plot_periodogram(x,y,N_freq,min_freq,max_freq,height_threshold,"periodogram")
plot_periodogram(x,y,N_freq,min_freq,max_freq,height_threshold,"ls")

内容的提问来源于stack exchange，提问作者Prav001