question: What is the relation between Bockstein homomorphism and Steenrod square?
For example, can one explain why the following relation works in the case of cohomology group with $\mathbb{Z}_2$ coefficient?
For $x \in H^m(M^d,\mathbb{Z}_2)$, $$ \beta_2 x=\frac{1}{2} d x \text{ mod } 2 \in H^{m+1}(M^d,\mathbb{Z}_2) $$ which is the Bockstein homomorphism. It turns out that the Bockstein homomorphism is the same as the Steenrod square $$ \text{Sq}^1(x) =\beta_2 x. $$
Why does the relation holds?
This is proven in Hatcher's book as Theorem 4L.12. It would be way too long to give a full proof here, because I would basically need to redefine the Steenrod squares. But the basic idea is as follows.
Consider the generator $\omega \in H^1(\mathbb{RP}^2; \mathbb{Z}/2\mathbb{Z})$. One of the properties of $\operatorname{Sq}^1$ says that $\operatorname{Sq}^1(\omega) = \omega^2$, and you can also explicitly compute $\beta(\omega) = \omega^2$. There's nothing else on which the two operations could be nonzero, so they agree on $H^*(\mathbb{RP}^2; \mathbb{Z}/2\mathbb{Z})$; by stability, they agree on the cohomology of iterated suspensions of $\mathbb{RP}^2$, i.e. on $H^*(\Sigma^n \mathbb{RP}^2; \mathbb{Z}/2\mathbb{Z})$.
But $\Sigma^n \mathbb{RP}^2$ is the $(n+2)$-skeleton of the Eilenberg--MacLane space $K(\mathbb{Z}/2\mathbb{Z}, n+1)$. Using Brown's representability theorem, any class $x \in H^{n+1}(X; \mathbb{Z}/2\mathbb{Z})$ is the pullback of some class $$y \in H^{n+1}(K(\mathbb{Z}/2\mathbb{Z}, n+1); \mathbb{Z}/2\mathbb{Z}) \cong H^{n+1}(\Sigma^n \mathbb{RP}^2; \mathbb{Z}/2\mathbb{Z}).$$
Since both $\operatorname{Sq}^1$ and $\beta$ commute with pullbacks and agree on $y$, you get $\operatorname{Sq}^1(x) = \beta(x)$.