Unit impulse response $h[n]$ of a discrete time system with multiple characteristic roots 0

I am trying to obtain the unit impulse response of a system in the form: $$y[n+N]+...+a_{N-1}y[n+1]+a_{N}y[n]=b_0x[n+N]+...+b_{N-1}x[n+1]+b_{N}x[n]$$ $$Q[E]y[n]=P[E]x[n]$$ (where E is the unit advance operator) Using the following formula: $$h[n]=\frac{b_N}{a_N}\delta[n]+y_c[n]u[n]$$ Or if $$a_N=0$$ $$h[n]=A_0\delta[n]+A_1\delta[n-1]+y_c[n]u[n]$$ Where $$A_0$$ and $$A_1$$ are determined through iterative calculation of initial values of $$h[n]$$ and $$y_c[n]$$ is a combination of the systems characteristic roots that satisfies $$Q[E]y_c[n]u[n]=0$$

The system in question is: $$y[n]=0.5y[n-1]+x[n]-(0.5)^4x[n-4]$$

Upon converting the system to using the operator E notation to obtain the characteristic roots, I obtain 4 roots 3 of which are repeated zeros. $$(E^4-0.5E^3)y[n]=(E^4-0.5^4)x[n]$$ $$\gamma_{1,2,3} = 0, \gamma_4=0.5$$

Using the relations described above and by simply discarding the repeated zero roots I then obtain the answer: $$h[n] = (0.5)^nu[n]$$ However, when replacing $$x[n]$$ and $$y[n]$$ with $$\delta[n]$$ and $$h[n]$$ respectively in the original system difference equation, I obtain: $$h[n] = (0.5)^n(u[n]-u[n-4])$$ Which is matching with my model answer.

I figured that maybe there is a mistake in my original solution due to the discarding of the repeated zero roots. How should I approach solving systems with zero repeated roots?

If you have a difference equation

$$\sum_{k=0}^Na_ky[n-k]=\sum_{k=0}^Mb_kx[n-k]\tag{1}$$

and if $$M>N$$, the general expression for the corresponding causal impulse response is

$$h[n]=A_0\delta[n]+A_1\delta[n-1]+\ldots+A_{M-N}\delta[n-(M-N)]+u[n]\sum_{k}c_k\gamma_k^n\tag{2}$$

where $$\gamma_k$$ are the roots of the characteristic equation, and where I've assumed simple roots. Note that roots $$\gamma_k=0$$ do not contribute to $$(2)$$ (even if there are multiple roots equal to zero).

For the given example we have $$M-N=3$$ and one non-zero root $$\gamma=\frac12$$. Hence the solution has the form

$$h[n]=A_0\delta[n]+A_1\delta[n-1]+A_2\delta[n-2]+A_3\delta[n-3]+c\left(\frac12\right)^n\tag{3}$$

The constants $$A_i$$, $$i=0,1,\ldots,3$$, and $$c$$ must be determined by iterating through the difference equation with $$x[n]=\delta[n]$$, $$h[n]=y[n]$$, and $$h[-1]=0$$. If you do that you'll find that $$c=0$$, which means that the given difference equation describes an FIR system. In terms of transfer functions this means that there is a pole-zero cancellation, leaving only poles at the origin of the complex plane.