We consider the problem of finding a solution to a multivariate polynomial equation system of degree $d$ in $n$ variables over $\mathbb{F}_2$. For $d=2$, the best-known algorithm for the problem is by Bardet et al. [J. Complexity, 2013] and was shown to run in time $O(2^{0.792n})$ under assumptions that were experimentally found to hold for random equation systems. The best-known worst-case algorithm for the problem is due to Björklund et al. [ICALP'19]. It runs in time $O(2^{0.804n})$ for $d = 2$ and $O(2^{(1 - 1/(2.7d))n})$ for $d > 2$. In this paper, we devise a worst-case algorithm that improves the one by Björklund et al. It runs in time $O(2^{0.6943n})$ (or $O(1.6181^n)$) for $d = 2$ and $O(2^{(1 - 1/(2d))n})$ for $d > 2$. Our algorithm thus outperforms all known worst-case algorithms, as well as ones analyzed for random equation systems. We also devise a second algorithm that outputs all solutions to a polynomial system and has similar complexity to the first (provided that the number of solutions is not too large). A central idea in the work of Björklund et al. was to reduce the problem of finding a solution to a polynomial system over $\mathbb{F}_2$ to the problem of counting the parity of all solutions. A parity-counting instance was then reduced to many smaller parity-counting instances. Our main observation is that these smaller instances are related and can be solved more efficiently by a new algorithm to a problem which we call multiple parity-counting.