How can we calculate the gradient of the Boltzmann policy over reward function?
The main point here is that you can write $Q(s, a|\theta) = R = \theta^\top \phi(s, a)$. For more details on this, you can read up on Policy Gradients (Chapter 13) from the 2nd edition of the Sutton and Barto book (in fact the expression you're looking for is equation 13.9)

For simplicity, I'm setting $\beta=1$, but you can always put it in once you get the idea. Therefore, the expression for $\pi(s,a,\theta)$ (I'm including $\theta$ here to make the dependence of $\pi$ explicit) is now:
$$
\pi(s,a,\theta) = \frac{exp(\theta^\top \phi(s, a))}{\sum_{a'} exp(\theta^\top \phi(s, a'))} = \frac{A}{B}
$$

Also, it is standard to compute $\nabla_\theta\log \pi(s,a,\theta)$ instead of $\nabla_\theta \pi(s,a,\theta)$ (although you totally can), so I'll go ahead and do that. To differentiate $\log \pi(s,a,\theta)$ with respect to the parameters, $\theta$, let's do some calculus:
$$
\nabla_\theta \log \pi(s,a,\theta) = \nabla_\theta \log \frac{A}{B} = \nabla_\theta \log A - \nabla_\theta \log B \\
= \nabla_\theta \left[ \theta^\top \phi(s,a) \right] - \nabla_\theta \left[ \log \sum_{a'} exp(\theta^\top \phi(s, a')) \right] \\
= \phi(s,a) - \frac{1}{\sum_{a'} exp(\theta^\top \phi(s, a'))} \sum_{b} exp(\theta^\top \phi(s, b)) \phi(s,b) \\
= \phi(s,a) - \sum_{b} \frac{exp(\theta^\top \phi(s, b))}{\sum_{a'} exp(\theta^\top \phi(s, a'))} \phi(s,b) \\
= \phi(s,a) - \sum_{b} \pi(s,b,\theta) \phi(s,b)
$$

Therefore, you can write $\nabla_\theta\log \pi(s,a,\theta) = \phi(s,a) - \sum_{b} \pi(s,b,\theta) \phi(s,b)$