NSolve[x^5-16*x+1 ==0,x]

ord=5;
li=CoefficientList[Series[x^5-16*x+ε/.x->Sum[a[n]*ε^n,{n,0,ord}],{ε,0,ord}],ε];sol=Solve[li==0];r=Sum[a[n]ε^n,{n,0,ord}]/.sol

N[r/.ε->1]