What is the derivative of a Radial Basis Interpolation function? A radial basis interpolation function is described as:
$
f(\textbf{x})=\sum_{k=1}^N c_k \phi(\lVert \textbf{x}-\textbf{x}_k \rVert_2), \ \textbf{x}\in\mathbb{R}^s
$
where $\textbf{x}_k$ are the $N$ scattered points and $c_k$ are the coefficients of the function obtained by solving the linear system:
$
\begin{bmatrix}
f(\textbf{x}_1) \\ f(\textbf{x}_2) \\ \vdots \\ f(\textbf{x}_N) 
\end{bmatrix}
=
\begin{bmatrix}
\varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_1\rVert} - \frac{1}{\lVert\textbf{x}_{1}\rVert} \right\rvert\right)  &  \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_1\rVert} - \frac{1}{\lVert\textbf{x}_{2}\rVert} \right\rvert\right)  & \ldots & \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_1\rVert} - \frac{1}{\lVert\textbf{x}_{N}\rVert} \right\rvert\right)   \\ 
\varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_2\rVert} - \frac{1}{\lVert\textbf{x}_{1}\rVert} \right\rvert\right)  & \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_2\rVert} - \frac{1}{\lVert\textbf{x}_{2}\rVert} \right\rvert\right) & \ldots & \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_2\rVert} - \frac{1}{\lVert\textbf{x}_{N}\rVert} \right\rvert\right)  \\
\vdots & \vdots & \ddots & \vdots \\
\varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_N\rVert} - \frac{1}{\lVert\textbf{x}_{1}\rVert} \right\rvert\right)  & \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_N\rVert} - \frac{1}{\lVert\textbf{x}_{2}\rVert} \right\rvert\right)  & \ldots & \varphi\left(\left\lvert \frac{1}{\lVert\textbf{x}_N\rVert} - \frac{1}{\lVert\textbf{x}_{N}\rVert} \right\rvert\right)  \\
\end{bmatrix}
\begin{bmatrix}
c_1 \\ c2 \\ \vdots \\ c_N
\end{bmatrix}
$
The radial basis function $\phi$ can be a Gaussian, an inverse multiquadric, etc. For a Gaussian, we get:
$
\phi(x) = e^{-(\epsilon x)^2}
$
for some free parameter $\epsilon$
How do you find the derivative of the function and is it well-defined on all values or do we get into zero-denominator situations? for example, what is $\frac{df(x)}{dy}$ if $\textbf{x}\in\mathbb{R}^2$?

The coefficients $c_k$ are independent of $\mathbf x$. Therefore, the partial derivatives of $f$ are 
$$\frac{\partial}{\partial x_j}f(\textbf{x})=\sum_{k=1}^N c_k \frac{\partial}{\partial x_j}\phi(\lVert \textbf{x}-\textbf{x}_k \rVert_2)$$
By the chain rule, 
$$\frac{\partial}{\partial x_j}\phi(\lVert \textbf{x}-\textbf{x}_k \rVert_2) = 
\frac{\langle \mathbf e_j, \mathbf x-\mathbf x_k\rangle }{\|\mathbf x-\mathbf x_k\|_2}\phi'(\lVert \textbf{x}-\textbf{x}_k \rVert_2) 
$$
where $\mathbf e_j$ is the $j$th standard basis vector. That is, $\langle \mathbf e_j, \mathbf x-\mathbf x_k\rangle$ is just the $j$th component of $\mathbf x-\mathbf x_k$.
You can simplify this a bit by writing $\phi(x)=\psi(x^2)$. Then the argument of $\psi$ is squared norm, which is more pleasant to differentiate, the result being 
$$\frac{\partial}{\partial x_j}f(\textbf{x})=2 \sum_{k=1}^N c_k \langle \mathbf e_j, \mathbf x-\mathbf x_k\rangle \psi'(\lVert \textbf{x}-\textbf{x}_k \rVert_2)$$
In any event, $f$ has derivative everywhere, as long $f$ itself is defined. 
