张春成
2022/09/22阅读：30主题：默认主题
样本点影响范围的定量估计
样本点影响范围的定量估计
对于连续时间信号，我们总可以用离散时间的样本点来近似地逼近它。
由于离散时间采样和噪声的影响，需要采用复杂的计算方法，才能从样本点出发，对连续时间信号的原始分布进行估计。
由于我们假设信号是连续的，因此就不能把观测到的样本点孤立地看待，而是要认为这些样本点之间彼此具有某种关联。要分析这种关联，就需要对这些样本点的影响范围和强度进行估计。
本文是一个 DEMO，尝试通过 B 样条的方法对孤立样本点在连续时间轴上的影响范围和强度进行定量估计。
基本假设
首先假设有一个连续时间信号，我们对它进行采样。根据得到的采样值对原始信号进行估计，估计的准则有很多，在此采用 B 样条的估计方法。它在保持曲线连续性的基础上，使整体误差最小。
这些误差我们认为是离散时间采样和噪声造成的“系统误差”。也就是说，我们在这里臭不要脸的认为原始信号就是蓝色曲线的样子。
.png)
连续信号的采样及 B 样条还原
噪声模拟及分析
为了考察每个样本点对整体信号的“贡献度”，我们采用逐点加噪的方法，对采集到的离散时间样本点进行噪声处理。
目前增加噪声的方法比较简单粗暴，就是将它从观测值集合中去除。去除后，我们当然可以利用新的，不包含目标样本的子集对原始信号进行估计，估计结果如下。
.png)
估计值的点图
.png)
估计值的连线图
其中，红色（full）曲线是全部样本点估计得到的原始真实曲线，其他颜色的曲线分别是去掉不同样本得到的估计曲线。我们看到，单一样本的缺失，可能造成的影响大小不一。
这时，我们对于每一个样本点都能得到两条不同的曲线，
-
一条当是它存在时，估计出的原始信号曲线； -
另一条是当它不存在时，估计出的原始信号曲线；
从逻辑上讲，只要将两条曲线进行对应时间相减，得到的差值曲线就是当某个样本点不存在时，它的缺失对连续信号造成的影响大小。值越大说明该样本点在此处的影响力越强，反之则说明影响力越弱。差值曲线如下
.png)
差值曲线
这些曲线即是某个样本点在整条时间轴上的影响力曲线，可以看到它不仅与样本点所处的时间点位置有关，更加具有一定的“广延性”，即它的影响力大小同时受到它与周围样本点之间关系的影响。下图将这种影响画成类似协方差矩阵的视图形式，就可以看得更加明显一些
.png)
差值曲线的矩阵视图
从图中可以看到，有些样本造成的影响虽然在幅度上不是最高，但能够在时间上延伸很远的距离；并且这种延伸呈线准周期性的衰减模式；另一些样本点，虽然它的缺失所导致的局部估计值的变化幅度较大，但它的影响却局限在一小段时间内。也就是说，样本点的影响范围和幅度大小之间并不是简单的时间临近和线性衰减关系。
分析代码
import copy
import numpy as np
import pandas as pd
import plotly.express as px
import plotly.graph_objects as go
import geomdl
from geomdl import BSpline
from geomdl.visualization import VisPlotly
from geomdl.visualization import VisVTK
# -----------------------------------------
# Main setting
ctrlpts_num = 20
evaluate_num = 100
template = 'plotly_dark'
degree = 3
params = dict(
degree=degree,
delta=0.025
)
# -----------------------------------------
# Toolbox
def draw(fig, template=template):
fig.update_layout(template=template)
fig.show()
def estimate(params, render=False):
'''
Estimate the BSpline using params
'''
# Compute knotvector
n = len(params['ctrlpts'])
params['knotvector'] = geomdl.knotvector.generate(
params['degree'], n)
# Setup
curve = BSpline.Curve()
curve.degree = params['degree']
curve.ctrlpts = params['ctrlpts']
curve.knotvector = params['knotvector']
curve.delta = params['delta']
# Evaluate
curve.evaluate()
# Visualize
if render:
curve.vis = VisPlotly.VisCurve3D()
curve.render()
# Estimate
x = np.linspace(0.0, 1.0, evaluate_num)
df = pd.DataFrame(curve.evaluate_list(x), columns=['x', 'y'])
df['idx'] = range(evaluate_num)
df['group'] = params['group']
# Record
estimations[params['group']] = df
print('The estimations contains', [e for e in estimations])
return df
# -----------------------------------------
# Generate data
x = sorted(np.random.randn(ctrlpts_num, 1))
y = np.random.randn(ctrlpts_num, 1)
ctrlpts = np.concatenate((x, y), axis=1)
ctrlpts = ctrlpts.tolist()
print('Shape of ctrlpts', np.array(ctrlpts).shape)
# -----------------------------------------
# Prepare container
df = pd.DataFrame(ctrlpts, columns=['x', 'y'])
df['idx'] = range(ctrlpts_num)
df['group'] = 'ctrlpts'
estimations = dict(
ctrlpts = df
)
# -----------------------------------------
# Estimate - full
params['group'] = 'full'
params['ctrlpts'] = ctrlpts
df = estimate(params, render=True)
df
# -----------------------------------------
# Estimate - noise
for j in range(1, ctrlpts_num-1):
params['group'] = 'pop-{}'.format(j)
_ctrlpts = copy.copy(ctrlpts)
_ctrlpts.pop(j)
params['ctrlpts'] = _ctrlpts
df = estimate(params, render=False)
# -----------------------------------------
# Concat container
table = pd.concat([it for it in estimations.values()])
table
# -----------------------------------------
# Display
# Estimation
fig1 = px.scatter(table.query('group=="ctrlpts"'), x='x', y='y', title='Ctrl pts')
fig1.data[0]['marker']['color'] = '#eeeeee'
fig2 = px.line(table.query('group=="full"'), x='x', y='y', title='BSpline')
fig = go.Figure()
fig.add_traces((fig1.data[0], fig2.data[0]))
fig.update_layout(title='Estimation')
draw(fig)
# Estimation noise
df = table.copy()
df['size'] = 1
fig = px.scatter(df, x='x', y='y', color='group', title='Compare (scatter)')
draw(fig)
fig = px.line(df, x='x', y='y', color='group', title='Compare (line)')
draw(fig)
# -----------------------------------------
# Compare
# Compare in line
baseline = np.array(table.query('group == "full"')[['y']])
x = np.array(table.query('group == "full"')[['x']])
dfs = dict()
for p in range(1, ctrlpts_num-1):
group = 'pop-{}'.format(p)
d = np.array(table.query('group == "{}"'.format(group))[['y']])
diff = np.abs(d - baseline)
df = pd.DataFrame(diff, columns=['diff'])
df['x'] = x
df['group'] = group
dfs[group] = df
df = pd.concat([e for e in dfs.values()])
fig = px.line(df, x='x', y='diff', color='group', title='Effect of absent (Curve)')
draw(fig)
# Compare in mat
mat = np.array([np.array(e['diff']) for e in dfs.values()])
fig = px.imshow(mat, title='Effect of absent (Matrix)')
draw(fig)
作者介绍