Python绘制3×3互相关子图时遇ValueError: a cannot be empty的解决方法

问题描述 #

需要将互相关（cross-correlation）结果绘制为3×3的子图，但运行代码时出现错误「ValueError: a cannot be empty」，尝试多种方法无法解决，求修复方案。

原始代码 #

import numpy as np
import pandas as pd
np.random.seed(1)

rows,cols = 540,2
data = np.random.rand(rows,cols) 
tidx = pd.date_range('1975-01-01', periods=rows, freq='MS') 
tim=pd.to_datetime(tidx)
df_t=pd.DataFrame(tim, columns=['tim']) # time series 
df_d = pd.DataFrame(data, columns=['a','b']) # data set 
df_f=pd.concat([df_t, df_d], axis=1) # final dataframe 
t1_a=df_f['tim'].iloc[0]
N=np.linspace(0, 540, 10)
M=np.arange(9)
fig, axes = plt.subplots(3,3, figsize=(12, 10), sharex=True, sharey=True,squeeze=True)
plt.subplots_adjust(wspace=0,hspace=0)

T=('[a:1975-1980]', '[b:1980-1985]', '[c:1985-1990]', '[d:1990-1995]', '[e:1995-2000]', '[f:2000-2005]', '[g:2005-2010]', '[h:2010-2015]', '[i:2015-2020]')

for i, j, k in zip(N, M, T):
    
    t2_a=t1_a + relativedelta(months=i)
    t3_a=t1_a + relativedelta(months=i+60)
    sub=df_f[(df_f.tim>=t2_a) & (df_f.tim<t3_a)]
    t1_b=sub['tim'].iloc[0]
    ax=axes.flatten()[j]
    lag=[]
    cr=[]
    for m in range(0, 60, 6):
        t2_b=t1_b + relativedelta(months=m)
        t3_b=t1_b + relativedelta(months=m+12)
        sub_b=sub[(sub.tim>=t2_b) & (sub.tim<t3_b)]
        w_sub=sub_b[['a', 'b']]
        e_sub=sub_b[['a', 'b']]
        wl=w_sub['a'] - np.mean(w_sub['a'])
        eq=e_sub['b']- np.mean(e_sub['b'])
        wl=wl/np.max(wl)
        eq=eq/np.max(eq)
        corr= scipy.signal.correlate(wl-np.mean(wl), eq-np.mean(eq), mode=str(mode), method=str(method))
        lags = signal.correlation_lags(len(wl), len(eq), mode=str(mode))
        cr.append(corr)
        lag.append(lag)
        
    cc=pd.DataFrame(cr)
    ccc=np.mean(cc)
    lll=pd.DataFrame(lag)
    lll=np.mean(lll)
    ax.plot(lll, ccc, linewidth=1.0, color='black') 
    ax.text(0, 0.6, k, fontsize=14)
    fig.tight_layout()

#fig.savefig('Example.jpg', dpi=500)    

错误原因分析 #

1. 缺失必要导入：代码使用了matplotlib.pyplot、dateutil.relativedelta、scipy.signal但未导入，会引发基础运行错误。
2. 时间切片超出数据范围：N=np.linspace(0,540,10)包含540，此时i+60=600，超出数据集的540个月范围，导致sub为空；内层循环m到54时，m+12=66超出sub的60个月范围，导致sub_b为空。
3. 列表追加逻辑错误：lag.append(lag)会导致列表嵌套自身，最终计算均值时数据结构异常。
4. 未定义参数：mode和method变量未赋值，scipy.signal.correlate调用会报错。
5. 除以零风险：若wl或eq的最大值为0，执行wl/np.max(wl)会触发除以零错误。

修复方案及完整代码 #

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from dateutil.relativedelta import relativedelta
import scipy.signal as signal

np.random.seed(1)

rows, cols = 540, 2
data = np.random.rand(rows, cols)
tidx = pd.date_range('1975-01-01', periods=rows, freq='MS')
tim = pd.to_datetime(tidx)
df_t = pd.DataFrame(tim, columns=['tim'])  # 时间序列
df_d = pd.DataFrame(data, columns=['a', 'b'])  # 数据集
df_f = pd.concat([df_t, df_d], axis=1)  # 合并后的DataFrame
t1_a = df_f['tim'].iloc[0]

# 修正N的取值，确保i+60 <= 540，避免切片超出范围
N = np.linspace(0, 540-60, 9)
M = np.arange(9)
fig, axes = plt.subplots(3, 3, figsize=(12, 10), sharex=True, sharey=True, squeeze=True)
plt.subplots_adjust(wspace=0, hspace=0)

T = ('[a:1975-1980]', '[b:1980-1985]', '[c:1985-1990]', '[d:1990-1995]', '[e:1995-2000]', '[f:2000-2005]', '[g:2005-2010]', '[h:2010-2015]', '[i:2015-2020]')
# 定义互相关的mode和method参数
mode = 'full'
method = 'auto'

for i, j, k in zip(N, M, T):
    t2_a = t1_a + relativedelta(months=int(i))
    t3_a = t1_a + relativedelta(months=int(i)+60)
    sub = df_f[(df_f.tim >= t2_a) & (df_f.tim < t3_a)]
    # 检查sub是否为空，避免后续报错
    if sub.empty:
        print(f"子数据集为空，跳过索引{j}")
        continue
    t1_b = sub['tim'].iloc[0]
    ax = axes.flatten()[j]
    lag_list = []
    cr_list = []
    
    for m in range(0, 60, 6):
        t2_b = t1_b + relativedelta(months=m)
        t3_b = t1_b + relativedelta(months=m+12)
        sub_b = sub[(sub.tim >= t2_b) & (sub.tim < t3_b)]
        # 检查sub_b是否为空
        if sub_b.empty:
            print(f"子时间段{m}-{m+12}无数据，跳过")
            continue
        
        w_sub = sub_b['a']
        e_sub = sub_b['b']
        
        # 计算去均值后的序列
        wl = w_sub - np.mean(w_sub)
        eq = e_sub - np.mean(e_sub)
        
        # 避免除以零，若最大值为0则跳过归一化
        wl_max = np.max(np.abs(wl))
        eq_max = np.max(np.abs(eq))
        if wl_max != 0:
            wl = wl / wl_max
        if eq_max != 0:
            eq = eq / eq_max
        
        # 计算互相关和滞后
        corr = signal.correlate(wl - np.mean(wl), eq - np.mean(eq), mode=mode, method=method)
        lags = signal.correlation_lags(len(wl), len(eq), mode=mode)
        
        cr_list.append(corr)
        lag_list.append(lags)
    
    # 确保有有效数据再绘制
    if cr_list and lag_list:
        cc = pd.DataFrame(cr_list)
        ccc = np.mean(cc, axis=0)
        lll = pd.DataFrame(lag_list)
        lll = np.mean(lll, axis=0)
        ax.plot(lll, ccc, linewidth=1.0, color='black')
    
    ax.text(0, 0.6, k, fontsize=14)
    fig.tight_layout()

# fig.savefig('Example.jpg', dpi=500)
plt.show()

关键修复点说明 #

• 补全所有必要的导入语句，确保代码能正常运行。
• 调整N的取值范围，保证每个时间切片都能获取到60个月的有效数据。
• 添加空数据集检查，避免对空数据执行计算操作。
• 修正lag列表的追加逻辑，改为lag_list.append(lags)，确保数据结构正确。
• 定义mode和method参数，明确互相关的计算模式。
• 增加除以零的判断，避免归一化步骤出错。

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