How are electrons really moving in an atom? Niels Bohr proposed the solar system model of the atom, which is the most recognizable, and assumed electrons revolve in circular paths called orbits/energy levels. However, we know that the Bohr model is incorrect due to later observations, like the Stark and Zeeman effect.
Now the quantum model of an atom is the best model which we have so far which tells us that electrons are located in regions of probability called orbitals.

Are electrons really in motion or what is the nature of their motion inside the atom? Is the motion of an electron predictable or observable by some means or is it really haphazard motion, which is unpredictable?
EDIT:
Thank you everyone for kind attention, all answers really helped me understand and clear my confusion.
 A: While classical mechanics describes the motion of objects, quantum mechanics  describes the object in terms of a probability amplitude density. Thus, quantum mechanics assumes that the object takes all possible paths, which are consistent with the information available at the start of the experiment. Therefore, the motion of an electron inside an atom is not described by our modern theory.
A: There is no such thing as a "trajectory" of a particle in quantum mechanics. If a particle followed a single trajectory -- that is, for every time $t$, you could associate a unique position $x(t)$ where it was located -- then you could also assign a unique momentum $p(t)=m\dot{x}(t)$. This would contradict the Heisenberg Uncertainty Principle, which says you cannot know the position and momentum to arbitrary precision simultaneously.
The orbitals you plotted represent the probability of finding an electron at a given location in the atom for states of the atom with definite energy. That is all the information that exists about the position of the atom; in quantum mechanics, there is no "extra information" about where the particle is beyond what you see in the wave function, until you make a measurement of position.
A: From your question I take that you have at least a basic of QM but if anything is unclear feel free to ask in the comments.
According to many interpretations of QM, a particle is the same as its wavefunction. So if a wavefunction is not in motion, the particle isn't either. Although there are some caveats to this.
The wavefunction itself can move though As you can see in this animation of a 1D wavepacket inside an infinite square well. It shows $\psi(x)$ with the x-axis from left to right and the other two axes are the real and imaginary parts of the wavefuntion.

One way to describe this movement is by looking at points that have the same phase. In the animation above you can track such a point in your imagination and see that it follows the overall motion of the wavepacket. I will come back to this but first I have to consider something that is important for discussing orbitals: the difference between real and complex atomic orbitals (see also here).
In the absence of a potential the stationary solutions to the Schrödinger equation are plane waves of the form
$$\psi_k(x)=A_ke^{ikx}+B_ke^{-ikx}.$$
This forms a basis and we can write a general solution as a superposition of these basis elements. Using superposition we can also write this solution as follows
$$\psi_k(x)=C_k\sin(kx)+D_k\cos(kx).$$
I will call this second basis the 'real' basis for obvious reasons. Note that each coefficient can generally still be complex.
For the atomic orbitals there is a similar split between 'real' orbitals and 'complex' orbitals. The real orbitals are often used for chemistry and the complex orbitals more often in physics. The picture you show shows real orbitals. Using this link by Falstad we can animate the complex orbitals. In the animations below $n=1,l=1$ and $m=0$ (left) and $m=1$ (right). Note that the left picture is shown from the side while the right one is top-down. From the right animation it finally becomes clear what quantum angular momentum has to do with rotations. The wavefunction itself is rotating! The probability density $|\psi|^2$ is still stationary, only the phase is rotating. This means that if an electron is in this state then we would call this electron stationary. The average position is constant $\langle r\rangle=\text{constant}$.


A: 
So my doubt is that are electrons really in motion or what is nature of their motion inside the atom

Electrons are quantum mechanical entities. Quantum mechanics was invented in order to mathematically model the spectra  of atoms, black body radiation, the photoelectric effect to start with, which could not be mathematically modeled with the classical theories.
We know electrons as particles obeying classical equations when they are not bound in a potential, and their interaction with matter being  so small as to show a track that can be calculated with classical equations, as this electron in a bubble chamber was hit off a proton and its track identifies it as a "classical" particle, the small bubbles of its interaction with the medium tracing its path that can be fit with classical electrodynamic equations.

But an electron bound in an atom cannot be described with those equations, the basic reason being that orbiting about the nucleus it would classically radiate and lose energy and fall on the nucleus neutralizing it. That is why the Bohr atom has the axiom of angular momentum quantization in order to explain the stability of atoms.
Quantum mechanics with its wavefunction and its postulates is a theory, that fits the data up to now very well. With this theory only the probability of finding an electron at (x,y,z,t) can be calculated, connected with the wave function.

Is the motion of electron predictable or observable by some means or is it really haphazard motion which is unpredictable?

Within the atom only the probability distributions are predictable and testable by experiment. So it is not haphazard/random  but weighted by the $Ψ^*Ψ$ '
A: The quantum mechanical model of an atom is a fantastic calculation tool. The  electron is modeled as a point particle in a superposition of position eigenstates that adds to a wave function that can be an eigenstate of energy, angular momentum and a projection of angular momentum on an axis (ignoring spin).
If the full theory, the electron is a quanta of the electron field. There is no need to consider it a particle with a probabilistic position at all times, it's just that in quantum mechanics we consider it a point particle and that makes the description of an orbital challenging using classical intuition. The idea that it is moving around is not useful....unless quantum numbers become large and the wave function starts to approach a classical orbit.
With all those caveats out of the way: an eigenstate of the Coulomb potential is a stationary state: the time evolution is a global phase factor, which is not a physical observable. That is: nothing ever changes, there is no motion. An $l=1$ state is forever stuck in $l=1$ (perturbations do change this an allow transitions).
The closest thing to velocity/motion is the probability current density:
$${\bf j} =\frac 1 {2m}\big(\psi^*\hat{\bf p}\psi -  \psi\hat{\bf p}\psi^*\big)$$
which will show probability orbiting the nucleus in $l\ne 0$ states, but in terms of $\psi({\bf r})e^{iE_nt/\hbar}$, nothing about $ \psi({\bf r})$ is changing.
A: 
Are electrons really in motion or what is the nature of their motion
inside the atom?

If we take the correspondence principle seriously, we probably have to conclude that electrons are in motion inside the atom.

Is the motion of an electron predictable or observable by some means
or is it really haphazard motion, which is unpredictable?

I would say the notion of predictability is ambiguous.
Is classical mechanics predictable? In principle, yes, if one can determine precise initial coordinates and velocities and integrate the equations of motions. In practice, sometimes it is possible, sometimes it is not. And some classical systems are chaotic.
In classical electrodynamics, one can make approximate predictions, but one cannot know the initial electromagnetic field in every point.
In quantum mechanics, one has equations of motion (say, for the wave function of one or several particles) and can integrate them in principle, but it is problematic to determine the initial wave function. Furthermore, if the relevant classical system, such as that in the three-body problem, is chaotic, the quantum system, such as the helium atom, is chaotic as well (see Scholarpedia, 8(4):9818, where the spectum of the helium atom is calculated based on classical trajectories. By the way, this relatively recent development shows that the Bohr model can be expanded to more complex systems than the hydrogen atom).
On a different note, let me mention my own work (Entropy 2022, 24(2), 261; Quantum Rep. 2022, 4(4), 486-508), where I show that electron orbitals can be approximated arbitrarily well by plasma-like collections of electrons and positrons, and each of these particles has a trajectory, but the uncertainty principle still holds true for the collection as a whole.
A: The answer to this question depends to some degree on your definition of 'motion'. There is a 'wave' around the nucleus, the values of which are complex numbers. For the standard basis of angular momentum eigenstates, the magnitude of these complex numbers is fixed over time, but the phases for many of the states rotate around the nucleus. Since global changes of phase make absolutely no difference to the physics (phase differences do have an effect, of course), it is questionable whether this constitutes anything physical 'moving', but if we treat the wavefunction phase as something 'real', then yes, you get waves spinning around the nucleus, just as Bohr thought.
The usual treatment separates out the time-dependence by solving the time-independent Schrodinger equation. The actual time-dependent wavefunction $\psi(t)$ is a complex exponenential of time multiplied by a fixed time-independent function $\psi_0$.
$$\psi(t,r,\theta,\phi) = e^{-iEt/\hbar}\psi_0(r,\theta,\phi)$$
And the $\psi_0$ has a complex exponential dependence on azimuth, $\phi$. For example, consider the $2p$-state:
$$\psi_{2,1,\pm 1}=\mp\frac{1}{8\sqrt{\pi}a_0^{3/2}}\frac{r}{a_0}e^{-r/2a_0}\sin\theta\, e^{\pm i\phi}$$
The final $e^{\pm i\phi}$ combines with the $e^{-iEt/\hbar}$ to yield a rotating movement of the phase in azimuth $e^{i(-Et/\hbar \pm \phi)}$, the rate of rotation being related to the energy of the orbital.
For higher order orbitals, the spherical harmonic part of the wavefunction is of the form:
$$Y^m_l=N\,P^m_l(\cos\theta)e^{-im\phi}$$
where $P^m_l$ is the associated Legendre functions which gives the dependence on polar angle (colatitude), and $e^{-im\phi}$ is the azimuthal dependence which, for order $m\neq 0$ rotates around the complex unit circle as discussed. (For $m=0$, the wavefunction sort of sloshes backwards and forwards along an axis; like a planet that falls directly towards a star and out the other side, then back again.)
This all works in exactly the same way as for pure basis states of linear momentum, for which the intuition may be better developed. A momentum eigenstate is a complex exponential for which the phase values move through space at a fixed velocity, called the phase velocity. However, the complex magnitude is constant, everywhere, and so the particle has no determined position. It is everywhere at once, equally, and continues to be so for all time. It is completely unlocalised, and so nothing appears to change as it moves.
We only get observable 'motion' when we combine momentum eigenstates to create a localised wavepacket, the waves of a range of wavelengths combining to yield a moving lump of probability. This lump moves at a different velocity to the phases, called the group velocity.
In just the same way, we can form angular momentum wavepackets that observably rotate around the nucleus by superposing pure states. For example, the posts here and here show how the superposition of two eigenstates yields an oscillating electron probability.
The pure angular momentum eigenstates behave very much like the pure linear momentum eigenstates, in that only the unobservable complex phase 'moves'. Complex magnitudes and hence probabilities are constant over all azimuths/positions, so motion makes no difference. And combining basis eigenstates to create a lump of probability results in motion at a completely different velocity to any of its constituents. Nevertheless, to whatever extent you consider that the complex phase is 'physically real', there is a real circulating motion around the nucleus associated with angular momentum eigenstates.
Those orbital diagrams are possibly a little misleading if they're not explained properly. I think they're generally doing something like plotting the real part of the complex wavefunction, so you can see the periodicity and structure. But the azimuthal dependence being of the form $e^{\pm im\phi}$ means that we are only seeing part of the full picture.
