What are the rules behind vector product in gradient?
It helps to think of each output dimension separately.
You have $X$ which is an $(2 \times 1)$ vector, and $W_1$ is a $(1 \times 2 )$ vector. Their product is a scalar, of which we then take the sigmoid to get our output $Y_1$.
The gradient of this w.r.t. $W_1$ will be $f(WX) (1 - f(WX)) X^T$, which has the appropriate dimensions.
Then, you just stack these for all your outputs, giving you the shape you got.