I was reading the Mathworks documentation about the CSAPS and TPAPS. I got confused in 2D case where I can't understand how these two penalties are different? Hows the PSS looks like in the case of 2D for CSAPS?
csaps creates a tensor-product smoothing spline, which means that the roughness penalties in the two (or more) directions are handled independently of each other. So in the 2D case, the roughness penalty would have a term for $|D_1D_1f|^2$ and a term for $|D_2D_2f|^2$, each with their own smoothing parameters
p that can set independently of one another.
On the other hand,
tpaps (for 2D data only) creates a thin-plate smoothing spline, which means that the roughness penalty is based on the sum $|D_1D_1f|^2 + 2|D_1D_2f|^2 + |D_2D_2f|^2$. Note the additional term here for the mixed derivative. Also, in this case there is a single smoothing parameter
p that controls all three terms simultaneously.
An advantage of
tpaps is that its roughness penalty is invariant under rotations. This could be desirable if you are dealing with data where the two dimensions have the same units, such as if the function $f$ is supposed to represent a physical surface. In this case, there is a physical interpretation of the fit as (approximately) minimizing the surface's bending energy, subject to the constraint on the amount of error between the fit and the data. On the other hand,
csaps could be more appropriate if the two dimensions are on different scales.