向量化插值matlab

我有一組N個在S個時間點測量的測量值（對於不同的測量，時間點是不同的）。我有兩個矩陣：向量化插值matlab

N - 表示的測量值

t時NXS矩陣 - 表示測量時間

一個NXS矩陣

我想生成表示線性內插的矩陣VI時間測量TI。非量化的代碼版本如下：

tic; 
VI = zeros([size(V,1), size(TI,2)]); 
for j = 1:size(V,1) 
    VI(j,:) = interp1(T(j,:),V(j,:),TI); 
end 
toc;

我想改寫這個代碼，以消除for循環，使得它使用矩陣操作和功能來實現。它可以被矢量化嗎？

來源

2012-09-19 Marc

請提供一些數據以改進您的解決方案。該配置文件將使我們能夠了解瓶頸在哪裏。 –

@Andrey - 對不起。這個問題的重點不是「我如何加快我的代碼？」關鍵是 - 「這可以通過使用矩陣運算/函數來完成」來消除for循環。我將編輯該問題以刪除對執行時間的引用，這是一種紅鯡魚。 – Marc

你看，循環不是邪惡的。你說表演不是你的主要興趣。你想要矢量化的原因是什麼？你想爲讀者提供更清晰的代碼嗎？也許爲了節省內存空間？ –

我在工作，所以沒有時間熟悉我自己interp1（我從來沒有用過它），所以如果下面的答案是不恰當的，我提前道歉。

話雖如此，因爲你的循環的迭代不依賴於彼此，矢量化應該是可能的。在我看來，你可以使用mat2cell和cellfun擺脫顯式循環。

我的意思是如下一個簡單的例子：

NumRow = 4; 
NumCol = 3; 
V = randn(NumRow, NumCol); 
VCell = mat2cell(V, ones(NumRow, 1), NumCol); 
A = cellfun(@sum, VCell);

當然，我做了什麼，相當於sum(V, 2)。但我認爲我使用的方法可以擴展到您的情況。 mat2cell函數將V轉換爲單元列向量，其中每個單元格包含一行V.然後，對cellfun的調用將sum函數應用於VCell中的每個單元格，並返回A中的結果。您可以簡單地將@sum替換爲@interp1，當然，也可以適當地將輸入調整爲cellfun。

讓我知道如果你不能得到它的工作，我會嘗試和一些更明確的事情，一旦我回家。另外，如果你確實得到它的工作，但它不會加快速度，我會很感興趣知道這一點，所以請在這裏發佈結果。

來源

2012-09-20 03:02:56

確實使用cellfun代表了對for循環的改進嗎？ MATLAB是否聰明並且與cellfun並行？ – Marc

@Marc值得一提的先生！我不知道。這就是爲什麼我說我會有興趣知道如果使用它加快了事情:-)所以我很好奇，並做了一些基本的測試。結果很有趣。我將把它們作爲一個新的SO問題發佈。 –

@Marc我對這個主題新的SO問題是[這裏]（http://stackoverflow.com/questions/12522888/arrayfun-can-be-significantly-slower-than-an-explicit-loop-in-matlab-why ） –

Matlab的在列，而不是行包數據，所以我懷疑你會看到速度的提高只是從周圍的循環變化，不過去行去在列：

[N, S] = size(V); 

VT = V';        % Value series in columns 
TT = T';        % Time series in columns 
VIT = zeros(length(TI), N);   % Interpolated value series in columns 

for j = 1:N 
    VIT(:, j) = interp1(VT(:, j), TT(:, j), TI); 
end 

VI = VIT';       % Interpolated value series in rows

不幸的是我不知道由於不同的價值序列沒有相關的時間序列，因此認爲可以做更多的事情來進一步提高此代碼的性能。如果他們有同樣的時間，這樣我們就可以崩潰T與length(T) = S載體，那麼這將是很容易做到這一點：

VIT = interp1(VT, T, TI);   % (see VIT and VT from above)

來源

2012-09-20 04:45:37

試過這個;沒有加速執行;但感謝您的建議！ – Marc

如果你有，就像你說的，所有測量不同的測量時間（T是一個矩陣，而不是一個向量），你可以做你想做的一個調用arrayfun如下：

VI = arrayfun(@(x)(interp1(T(x,:),V(x,:),TI)), 1:size(V, 1), 'UniformOutput', false); 
VI = cell2mat(VI');

arrayfun類似於循環，但由於它是一個內部的MATLAB函數可能是更快。它返回一個矢量單元，所以第二行確保你有一個矩陣作爲輸出。你可能不需要它 - 這取決於你後來對數據做了什麼。如果另一方面對於不同的N值（T是尺寸爲S的向量，而不是矩陣，或者換句話說，T的所有行相等），在相同的時間進行測量，您可以插入一個調用interp1。

VI = interp1(T(1,:), V', TI)

在這裏你必須轉置V，因爲interp1在列內插值。這是因爲MATLAB按列存儲矩陣（列在內存中是連續的）。如果將V作爲SxN矩陣傳遞，它可能允許interp1的更高效的並行化，因爲所有的CPU都可以以更高效的方式訪問內存。因此，我建議你在整個代碼中轉置你的矩陣，除非你因爲性能原因在別的地方依賴這個確切的數據佈局。

編輯由於矩陣的列布局，您可以通過轉置矩陣和逐列工作來改進原始代碼。以下版本在我的計算機上N = 1000，S = 10000和TI爲10000個元素的速度提高了大約20％。由於更高效的內存帶寬利用率，它可能會隨系統大小而增長。

tic; 
VI = zeros(size(TI,2), size(V,2)); 
for j = 1:size(V,2) 
    VI(:,j) = interp1(T(:,j),V(:,j),TI); 
end 
toc;

來源

2012-09-20 07:39:07 angainor

嘗試了數組函數;與for循環相比，並沒有加快執行速度。感謝您的建議！ – Marc

@Marc你可以說速度是你的問題的目標。查看更新答案。 – angainor

這是很難說沒有數據和運行探查什麼，但是如果你的數據進行排序，您可以使用interp1q代替interp，不進行數據的任何檢查。

從MATLAB的幫助摘自：

對於interp1q正常工作，x必須是單調遞增列向量。 Y必須是長度爲（x）行的列向量或矩陣。 xi必須是列向量

來源

2012-09-20 07:50:15

感謝您的建議; interp1q取代15％的執行時間。 – Marc

向量化插值matlab

回答

相關問題