In matrix form, we're looking for coefficients [tex]b_3,b_2,b_1,b_0[/tex] such that
[tex]\underbrace{\begin{bmatrix}3^3&3^2&3^1&3^0\\1^3&1^2&1^1&1^0\\(-1)^3&(-1)^2&(-1)^1&(-1)^0\\6^3&6^2&6^1&6^0\\7^3&7^2&7^1&7^0\end{bmatrix}}_{\mathbf A}\underbrace{\begin{bmatrix}b_3\\b_2\\b_1\\b_0\end{bmatrix}}_{\mathbf x}=\underbrace{\begin{bmatrix}4\\2\\1\\5\\9\end{bmatrix}}_{\mathbf b}[/tex]
The best-fit solution is given by [tex]\mathbf x=(\mathbf A^\top\mathbf A)^{-1}\mathbf A^\top\mathbf b[/tex]. You should end up with
[tex]\mathbf x=\begin{bmatrix}b_3\\b_2\\b_1\\b_0\end{bmatrix}\approx\begin{bmatrix}0.0495\\-0.3463\\0.9157\\2.1149\end{bmatrix}[/tex]
Attached is a plot of the given points and the best-fit solution.
