Obtaining polynomial coefficients for a 3-D spline with data that are not a meshgrid

I have 200 of (y, x_1, x_2) and I'd like to fit a cubic smoothing spline of y on x_1 and x_2 to obtain the estimates of a_0, a_1, a_2 etc. where y=a_0+a_1*x_1+a_2*x_1*x_2 ... The csaps command does this for a 2-D spline, however for higher dimensions it requires input in the form of a meshgrid of (x_1, x_2) i.e. I need y values for every of meshgrid(x_1,x_2), which is not how my data come. The command scatteredInterpolant seems to take (y, x_1, x_2) as they come (not necessarily a meshgrid) but it does not split out coefficient estimates as output but instead just finer, interpolated points. Tried interp, interp2, spline etc as well but they don't work either. Is there a way to obtain coefficients of the polynomial of a 3-D spline in Matlab? Thank you!