基于二维高斯核估计bwt与age联合密度的代码报错求解

1. 输入查询点t形状不匹配：你当前要计算二维变量的联合密度，每个查询点必须是二维坐标，而你生成的t是长度为50的一维数组，形状为(50,)，无法和形状为(1173,2)的样本点x做广播运算，直接导致报错。
2. 二维高斯核公式错误：d维高斯核的归一化系数为(2π)^(-d/2)，二维场景下应该是1/(2π)，你写的是一维高斯核的系数1/√(2π)，计算结果会有偏差。
3. 查询点取值范围不合理：你选的两个变量取值范围远大于你设置的t范围1.0~5.0，计算(t-x)/h时会得到绝对值很大的负数，exp运算后结果趋近于0，所以你会看到返回值全零。
4. 列选择不符合需求：你需求是计算bwt（出生体重）和age的联合密度，但代码里选的是gestation和age，需要修正列名。

import numpy as np
import pandas as pd

# 读取数据集
babies_full = pd.read_csv("https://www2.helsinki.fi/sites/default/files/atoms/files/babies2.txt", sep='\t')

# 按需求选择bwt和age两列
babies_full1 = babies_full[['bwt', 'age']]
x = np.array(babies_full1, 'int')
n, d = x.shape  # n是样本量，d是维度=2
h = 5  # 带宽

# 修正后的二维高斯核
def k_2dgauss(x):
    # 输入x形状为 (n, 2)，每个元素是 (查询点-样本点)/h
    return np.exp(-np.sum(x**2, axis=1) / 2) / (2 * np.pi)

# 修正后的核密度计算函数，支持单个查询点的密度计算
def mv_kernel_density_single(t, x, h):
    # t是单个二维点，形状为(2,)
    return np.mean(k_2dgauss((t - x) / h)) / (h ** d)

# 生成二维查询网格，范围匹配变量实际取值
bwt_min, bwt_max = x[:, 0].min(), x[:, 0].max()
age_min, age_max = x[:, 1].min(), x[:, 1].max()
# 每个维度生成50个点
bwt_grid = np.linspace(bwt_min, bwt_max, 50)
age_grid = np.linspace(age_min, age_max, 50)
# 生成网格坐标
BWT, AGE = np.meshgrid(bwt_grid, age_grid)
# 初始化密度矩阵
density = np.zeros_like(BWT)

# 遍历每个网格点计算密度
for i in range(BWT.shape[0]):
    for j in range(BWT.shape[1]):
        t = np.array([BWT[i, j], AGE[i, j]])
        density[i, j] = mv_kernel_density_single(t, x, h)

# 输出结果示例：前5行5列的密度值
print(density[:5, :5])

• 如果你需要加速计算，可以利用numpy广播机制替换循环，把网格点变形为(N, 2)的数组后批量计算，无需遍历每个点。
• 如果你确实要计算gestation和age的联合密度，把列选择部分改回['gestation', 'age']即可。

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