Tương quan chéo chức năng không định kỳ với NumPy

Question

Tôi có hai tập dữ liệu mà tôi đang cố gắng tương quan chéo. Chúng trông giống với chức năng arctan, vì vậy tôi đã sử dụng nó làm mô hình để tìm hiểu cách xử lý tín hiệu của mình.

x = linspace(-15, 15, 2**13) 
f1 = arctan(x) 
f2 = arctan(x + 2) 

Các câu hỏi tôi cần phải trả lời được, bao nhiêu để tôi cần phải thay đổi tín hiệu màu xanh lá cây để làm cho nó (chủ yếu) chồng chéo với một màu xanh? Tôi nghĩ nó sẽ đơn giản như việc tìm kiếm tối đa trong hàm tương quan chéo của f1 và f2 và tôi đã thực hiện theo lời khuyên ở đây: How to correlate two time series with gaps and different time bases?. Đây là những gì tôi đã cố gắng

c = correlate(f1, f2, 'full') 
s = arange(1-2**13, 2**13) 
dx = 30/2**13 
shift = s[c.argmax()]*dx 

Tôi mong chờ shift để bằng nhiều hay ít chính xác 2, nhưng trên thực tế nó chỉ 0.234. Điều này không có ý nghĩa với tôi; Tôi đã tìm thấy tọa độ x của mức tối đa trong tương quan chéo, nên được tìm thấy khi có chồng chéo tối đa của hai tín hiệu.

Bất kỳ ý tưởng nào về cách tính số lượng này cho loại chức năng này?

EDIT: Tôi cần thêm, cho dữ liệu thực tế của tôi, tất cả các giá trị này sẽ được giữa zero và một

EDIT EDIT: Các chức năng sau đây là thực sự giống như dữ liệu thật của tôi:

x = linspace(-15, 15, 400) 
f1 = (arctan(-x) + pi/2)/pi 
f2 = (arctan(-x + 2) + pi/2)/pi 

Vì vậy, sử dụng công thức đưa ra ở đây: http://paulbourke.net/miscellaneous/correlate/ tôi có thể viết một hàm tương quan chéo rằng miếng dữ liệu để chèn sang bên trái và bên phải số không:

def xcorr(x, y); 
    mx = x.mean() 
    my = y.mean() 
    sx = x.std() 
    sy = y.std() 
    r = zeros(2*len(x)) 

    for d in range(-len(x), len(x)): 
     csum = 0 
     for i in range(0, len(x): 
      yindx = i - d 
      if i - d < 0: 
       yval = 1 
      elif i - d >= len(x): 
       yval = 0 
      else: 
       yval = y[yindx] 
      csum += (x[i] - mx) * (yval - my) 
     r[d + len(x)] = csum/(sx * sy) 
    return r 

Với chức năng này, bây giờ tôi có thể làm

c = xcorr(f1, f2) 
s = arange(-400, 400) 
dx = 30/400 
shift = s[c.argmax()] * dx 

nào đi ra đến 2,025, đó là gần như là bạn có thể nhận được tới 2 với độ chính xác này. Vì vậy, có vẻ như Jamie là chính xác, vấn đề là làm thế nào numpy correlate hiện padding của tín hiệu.

Vì vậy, rõ ràng chức năng xcorr của tôi thực sự chậm khi đứng. Câu hỏi bây giờ là, có cách nào tôi có thể làm cho NumPy làm một điều tương tự, hay tôi nên tiến hành viết thuật toán của tôi trong C bằng cách sử dụng ctypes?

Answer 1

Giá trị tối đa của tương quan chéo là sự dịch chuyển mà tại đó tổng sản phẩm của hai tín hiệu của bạn là tối đa. Người ta sẽ nghĩ rằng sự tương quan chéo của hai tín hiệu, một sự thay đổi thời gian của tín hiệu kia, sẽ là tối đa ở ca làm việc. Và mặc dù điều này đúng với các tín hiệu vô hạn, nhưng không có hai tín hiệu dài, do đó, với hai tín hiệu dài bằng nhau, bạn chỉ có một tổng kết trên 2**13 giá trị khác 0 cho sự dịch chuyển bằng 0 và giá trị cao hơn từ kết hợp tốt hơn hai hàm, được bù trừ nhưng có ít giá trị khác 0.

Nếu tín hiệu của bạn bằng 0 ở +/- vô cùng, điều này sẽ không quá tệ. Vì nó là, tôi đã không thể đưa ra bất kỳ giải pháp hợp lý sử dụng tương quan chéo.

Answer 2

Như mọi người đã chỉ ra, mối tương quan chéo bị nhiễu bởi phần đệm ngoài dữ liệu.

Mặc dù cảm giác như bạn đang vứt bỏ dữ liệu tốt, nhưng tốt hơn là chỉ cắt bộ dữ liệu để có thể thực hiện tương quan mà không có giả định (ít nhất là khi so sánh với dữ liệu thực tế tương ứng với dữ liệu đã tạo cho padding).

x = linspace(-15, 15, 4000) 
f1 = (arctan(-x) + pi/2)/pi 
f2 = (arctan(-x + 2) + pi/2)/pi 

L4 = int(len(f2)/8) 
sf2 = f2[L4:-L4] 

c = correlate(f1-mean(f1), sf2-mean(f1), 'same') 
print "peak correlation occurs at:", x[argmax(c)] # -2.02925731433 

plt.plot(x, c) 
plt.show() 

Ngoài ra, tuy nhiên, tôi không chắc chắn rằng xcorr là phương pháp tốt nhất ở đây. Cách tổng hợp khoảng cách giữa các giá trị trục y cho các thay đổi khác nhau và lấy giá trị tối thiểu, tránh xa tất cả các vấn đề về vị trí số 0, ...