I found (from here) that current density can be found with the formula: $$\vec j= \frac{dI}{ds} \vec a$$ where $s$ is the cross sectional area centred on the point we are considering and with a normal vector in the same direction as the direction of the current , $\vec a$ which is a unit vector.
But I have no idea how to use it? How do we write $I$ in such a way that we can differentiate it with respect to $s$? Any examples, would be helpful.
Take for example a conducting body with an arbitrary shape and volume $V$ connected to a battery that will cause a current flow in the body. You don't know the microscopic structure of the body; so in principle for any point $ \mathbf{r} \in V $ you can have a different value of $I$: $$ I=I(\mathbf{r})=I(x,y,z) $$ Now, take any surface $ \Sigma $ in $V$; for the points $ \mathbf{r} \in \Sigma $, your current is a function of the surface $ \Sigma $, which you can express in parametric form: $$ \Sigma = \{ (x',y', z') \in V : x'=x'(u,v), y'=y'(u,v), z'=z'(u,v), \text{for proper values of } u \text{ and } v \} $$ and thus: $$ I(\Sigma) = I(x'(u,v),y'(u,v),z'(u,v))$$ Now you can differentiate the expression of $I$ with respect to the parameters $(u,v)$. For each $\Sigma$ you have different parametric equations, and differentiate gives you the derivative of $I$ with respect to that particular $ \Sigma$.
That's how I would see differentiating a function with respect to a surface; but I think the formula you reported has a bit of notation abuse and makes much more sense if you use: $$ \mathbf{J}= \rho\mathbf{v} $$ where $\rho$ is the volumetric charge density of the body and $v$ is the speed of the charges flowing in the direction $\mathbf{a}$. With this expression of $\mathbf{J}$ is easier to calculate the derivative.
That is a limiting procedure pretty much similar to the definition of derivative. Fix the point $x$ where you want to compute the current density and consider a small disk $D$ centred at that point. There will be a certain current $I(D)$ flowing across this disk. As this disk shrinks on to the point, the ratio $$\frac{\text dI}{\text d S}\Big\vert_x:=\lim_{\sigma(D)\to0} \frac{I(D)}{\sigma(D)},$$ where $\sigma(D)$ is the surface of the disk, will represent the density current at the point $x$, modulo the velocity vector $\mathbf a$. Observe that there is a geometric factor involved in this which depends how the (normal to the) surface is oriented w.r.t. $\mathbf a$ at the point $x$, which turns out to be the cosine of the angle between these two vectors.
