The Bessel functions of the first kind were first defined by Daniel Bernoulli, and then generalised by Friedrich Bessel. They are defined simply as the solutions to the following second order differential equations

$$x^2 \frac{d^2 y}{ dx^2 } + x \frac{dy}{ dx} + (x^2 – \alpha^2) y = 0, \ \ \ \ \alpha \in \mathbb{C},$$

where $\mathbb{C}$ denotes the set of complex numbers. These equations are now known as Bessel equations. Their solutions are called Bessel functions of the first kind and usually denoted by $J_{\alpha}(x)$.

In the especial cases when $\alpha$ is an integer, i.e. when we have

$$x^2 \frac{d^2 y}{ dx^2 } + x \frac{dy}{ dx} + (x^2 – n^2) y = 0, \ \ \ \ n \in \mathbb{N},$$

the Bessel functions $J_n(x)$ are also referred as cylinder functions or the cylindrical harmonics because they appear in the solution to Laplace’s equation in cylindrical coordinates.

To solve the Bessel equation one can use the Frobenius method which leads us (in a non entirely trivial way that you can find Here) to the following expression

$$J_{\alpha}(x) = \sum_{m=0}^{\infty}\frac{(-1)^m}{m!\Gamma(m+\alpha+1)} \left(\frac{x}{2}\right)^{2m+\alpha} ,$$

where $\Gamma(z)$ is the gamma function. Now that we know that the solutions exist, let’s take a look at them!

For any integer or positive $\alpha$, Bessel functions of the first kind are finite at $x=0$. For negative non-integer $\alpha$, Bessel functions of the first kind diverge as $x$ approaches zero.

Bessel functions are also valid for complex arguments $x$, and an important special case is that of a purely imaginary argument. In such case, the solutions to the Bessel equation are called modified Bessel functions (or hyperbolic Bessel functions) of the first kind. So, they are defined in terms of the Bessel functions of the first kind as:

$$I_{\alpha}(x) = i^{-\alpha} J_{\alpha}(ix).$$

**Bonus Random Fact: **The probability density function of a non-central chi-square random variable $\chi^2(k, \lambda)$ with $k>0$, degrees of freedom and non-centrality parameter $\lambda>0,$ can be written in terms of modified Bessel functions of the first kind as follows

$$f(x; k, \lambda) = \frac{1}{2} e^{-\frac{ (x+\lambda) }{2} } \left( \frac{x}{\lambda} \right)^{\frac{k}{4} – \frac{1}{2} } I_{\frac{k}{2} -1 } (\sqrt{\lambda x}).$$

### Python Implementation

Bessel functions of the first kind, as well as their modified versions, are implemented in the scypy library. Here is a snippet showing how to plot these functions.

```
# Load the libraries that we need
import matplotlib.pyplot as plt
import numpy as np
import scipy.special as sp
# Define a function to plot the Besself functions of the first kind
def plot_bessel_functions(alphas, xrange, modified=False, my_title=None, return_fig = False, **fig_kw):
x = np.linspace(xrange[0], xrange[1], 500)
fig, ax = plt.subplots(1, 1, **fig_kw)
for v in alphas:
if modified:
ax.plot(x, sp.iv(v, x), lw=1.5)
my_legends = ['$I_{'+str(n)+'}(x)$' for n in alphas]
else:
ax.plot(x, sp.jv(v, x), lw=1.5)
my_legends = ['$J_{'+str(n)+'}(x)$' for n in alphas]
ax.legend(my_legends)
if my_title:
plot_title=my_title
elif modified:
plot_title='Modified Bessel Functions of 1st Order'
else:
plot_title='Bessel Functions of 1st Order'
plt.title(plot_title)
if return_fig:
return fig
else:
plt.show()
# Now we can plot!
plot_bessel_functions(alphas=[-5, -3, 0, 3, 5], xrange=[-10, 10])
plot_bessel_functions(alphas=[-5, -3, 0, 3, 5], xrange=[-10, 10], modified=True)
```