1
$\begingroup$

I wish to estimate the noise of given trajectory from a measurement (data which is provided).

I consider that the noise is random gaussian noise.

I add here examples of two trajectories, one of the original data, and the other smoothened.

The trajectories are the measured velocity of some particles.

The model which is used to calculate the trajectories is quite complicated, but in the simplest form it is a langevin equation.

The results I would like to obtain at the end to have an estimation of the amplitude of the noise (considering it is an uncorrelated random gaussian noise).

I add here a csv file here

https://fastupload.io/o2gAHhJcEWRcxKT

Note that there are some NaN values in the file, and the plot is done by normalizing each time series with its average.

This is the plot code for the smoothened data in matlab

    figure
hold on
plot(t,smooth(vf,25)/mean(vf,'omitnan'),'-b'); 
plot(t,smooth(vb,25)/mean(vb,'omitnan'),'-r');    
xlabel('time'); ylabel('mean normalized intensity');  
box on
legend('vf','vb')

enter image description here

enter image description here

$\endgroup$
0

2 Answers 2

2
$\begingroup$

After reading in the data (and "correcting" the NaN values with the previous valid number):

import csv
with open('Q83666vfvb.csv', newline='') as csvfile:
    rdr = csv.reader(csvfile)
    next(rdr)
    data = list(rdr)
   
idx = 0
vf = np.zeros(len(data))
vb = np.zeros(len(data))
for tup in data:
    vf[idx] = vf[idx-1] if tup[1] == "NaN" else tup[1]
    vb[idx] = vf[idx-1] if tup[2] == "NaN" else tup[2]
    idx = idx + 1

I'd first try using a simple DC blocking filter:

$$\displaystyle y[n] = x[n] - x[n-1] + \alpha\, y[n-1]$$

from scipy import signal

alpha = 0.9

vf_dc = np.zeros(len(vf))
vb_dc = np.zeros(len(vb))
vf_dc[0] = 0
vb_dc[0] = 0
for idx in np.arange(len(vf)-1) + 1:
    vf_dc[idx] = vf[idx] - vf[idx-1] + alpha*vf_dc[idx-1]
    vb_dc[idx] = vb[idx] - vb[idx-1] + alpha*vb_dc[idx-1]

and then just take the standard deviation of the vf_dc and vb_dc measurements.

Comparing the two results (blue) with their originals (orange):

vf

vb


The next thing I'd look at doing is removing the "noise" component from the original signal:

vf_corrected = vf - vf_dc
vb_corrected = vb - vb_dc

and then run a first order smoother over them:

$$\displaystyle y[n] = \alpha\, y[n-1] + (1-\alpha) x[n]$$

alpha = 0.9

vf_corrected_smooth = np.zeros(len(vf_corrected))
vb_corrected_smooth = np.zeros(len(vb_corrected))
vf_corrected_smooth[0] = vf_corrected[0]
vb_corrected_smooth[0] = vb_corrected[0]
for idx in np.arange(len(vf)-1) + 1:
    vf_corrected_smooth[idx] = alpha*vf_corrected_smooth[idx-1] + (1-alpha)*vf_corrected[idx]
    vb_corrected_smooth[idx] = alpha*vb_corrected_smooth[idx-1] + (1-alpha)*vb_corrected[idx]
    
plt.figure(1)
plt.plot(vf_corrected_smooth)
plt.plot(vf_corrected)

plt.figure(2)
plt.plot(vb_corrected_smooth)
plt.plot(vb_corrected)

Corrected and smoothed vf

Corrected and smoothed vb


Note that the above is completely disregarding the other information you have given: that this data is generated from a system described by the Langevin equation.

If playing with the two alpha values in what I've done above doesn't give satisfactory results, the next step will be to do what Filipe Pinto says in their answer and write down the differential equations governing your measurements and derive the appropriate Kalman filter (or extended Kalman filter if the system ends up being nonlinear).

Please give us more detailed information about the model your measurements follow and we can probably derive the appropriate discrete-time Kalman filter for your problem. There is some literature already that seems to apply the Kalman filter equations to Langevin systems. Perhaps have a look at that if the simple approach I'm suggesting doesn't work for you.

$\endgroup$
3
  • $\begingroup$ Thanks, what would be your second order approximation as a filter? in case there is also a drift and not just diffusion. $\endgroup$
    – jarhead
    Jul 3, 2022 at 14:04
  • $\begingroup$ @jarhead Added some commentary to try to address your question. $\endgroup$
    – Peter K.
    Jul 3, 2022 at 17:46
  • 1
    $\begingroup$ The model is quite complicated as it involves a set of differential equations, however, for approximations, first order, or third order normal forms are mostly sufficient. I guess this will work the same for the filter. $\endgroup$
    – jarhead
    Jul 7, 2022 at 14:05
1
$\begingroup$

A Kalman filter or a particle filter could be an option if your trajectories can somehow fit an ARMA model. A dual Kalman filter is also a possibility, one to track trajectories and another to track the filter. The differences between the filter estimates and the signal you have would be your noise data.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.