Edit: as suggested by @Slereah the following explanation of the path integral approach to QED misses gauge fixing and ghost terms to be correct. So, how can this explanation be made correct?
What follows is a (hopefully partially correct) description of path integral approach to QED.
As per my understanding from reading some papers and books, the value of functional integral $Z$ is given by-
$$Z=\int_{-\infty}^{+\infty}\mathcal{D}\psi\mathcal{D}A\exp\left[\frac{i}{\hbar}\int_{\infty}^{+\infty}\mathcal{L}(\psi,A,\partial_{\mu}\psi,\partial_{\mu}A)d^{4}x\right]$$
Here, $$\mathcal{D}\psi=\sum\limits _{i}\mathcal{D}(\psi_{i}\phi_{i}(x)),$$ where $\psi_{i}$ is Grassmann number and $\phi_{i}$ is the basis of four-component spinors. Also, $A_{\mu}$ is the electromagnetic four potential. Finally,
$$\mathcal{L}=\bar\psi(i\gamma^\mu D_\mu-m)\psi -\frac{1}{4}F_{\mu\nu}F^{\mu\nu}$$
is the QED Lagrangian for the interaction between the spin half field and the electromagnetic field.
Is my understanding of the topic till now correct?
Now, last bit is the physical part. To me, this equations is saying that given two configurations of the $\psi$ and $A$ fields separated by time $T$, we can find the probability amplitude that this state change will happen by doing the functional integral of both the fields, spread everywhere in space and time, with only one constraint that the initial and final configurations of the field at times, let us say, $t_{1}$ and $t_{2}$ ( such that $t_{1}-t_{2}=T$) should be constant over the integral. Hence, the past as well as the future will affect the present, but most likely the effect will cancel out.
So, basically, we know the initial and final field configurations, and we try to find Lagrangian of all the field configurations spread everywhere in space and time which have our initial and final field configurations in common.
We put and integrate the values of these Lagrangian in the imaginary exponential function, which will result in most of the extreme field cases to cancel out, and the most natural and shortest ones to remain. This will give the probability amplitude and its square will give the probability for the event to happen.
So, is my physical interpretation correct? If no, then please point out how it can corrected...
I understand that the problem is not yet quite solved, because the integrals we are talking about is spread till infinities, and hence difficult to carry out. Also, some of the integrals will be on Grassmann variables, which are quite non intuitive. But do I get the basic principles correctly?