Ultrasound speckle tracking provides a robust motion estimation of fine tissue displacements along the beam direction. Extensions to 2-D have been proposed in recent years. Due to relatively coarse lateral sampling, several solutions relied on lateral interpolation in order to achieve subsample accuracy. We introduce a new multi-dimensional speckle tracking method (MDST) with subsample accuracy in all dimensions. The proposed algorithm is based on solving a least squares problem to estimate the coefficients of a second order polynomial expansion to fit the magnitude of the two dimensional complex normalized correlation of the generalized analytic signal in the vicinity of the true peak. The generalization method utilizes the Riesz transform which is the multidimensional Hilbert transform. The displacement is estimated from acquired successive radio-frequency data frames of the region of interest. Field II simulation of flow data in a channel with a bench mark known displacement is generated to validate the accuracy of the method. In addition, the new MDST method is applied to imaging data from a flow phantom (ATS Model 524) to estimate the flow motion and pulsating channel wall. Simulations and experimental results demonstrate the effectiveness of the proposed technique.