(another comment, not an answer.)

Experimentally, with the following code

from sympy import *
var('x')
var('n', integer = True)

f = ( (2+(x+3) * sqrt(-x))/(2*(x+4)) *
      sqrt(-x)**n + ((2-(x+3) * sqrt(-x)) / (2*(x+4))) * (-sqrt(-x))**n +
      + (x+2+sqrt(x*(x+4)))/(2*(x+4)) *
        ( (x+sqrt(x*(x+4)))/2 )**n +
        (x+2-sqrt(x*(x+4)))/(2*(x+4)) *
        ( (x-sqrt(x*(x+4)))/2 )**n
        )
pprint(f)
for i in range(10):
   eq = simplify(f.subs(n,i))
   print ('========== n = ', i)
   pprint(eq)
   print ('have solutions')
   sols = solve(eq)
   pprint (sols)
   pprint ('approx. = ')
   pprint ([s.evalf() for s in sols])

we see the following: $$ \begin{array}{l} f_{0} = 1\\ f_{1} = 0\\ f_{2} = x^{2}\\ f_{3} = x^{2} \left(x + 2\right)\\ f_{4} = x^{2} \left(x^{2} + 2 x + 1\right)\\ f_{5} = x^{3} \left(x^{2} + 3 x + 1\right)\\ f_{6} = x^{4} \left(x^{2} + 4 x + 4\right)\\ f_{7} = x^{4} \left(x^{3} + 5 x^{2} + 7 x + 3\right)\\ f_{8} = x^{4} \left(x^{4} + 6 x^{3} + 11 x^{2} + 6 x + 1\right)\\ f_{9} = x^{5} \left(x^{4} + 7 x^{3} + 16 x^{2} + 13 x + 2\right)\\ \end{array} $$

And it seems that already $f_9$ have some imaginary roots, albeit very small. Do these imaginary roots really converges to 0 when $n\to \infty$ or they stay on the same magnitude ?

(another comment, not an answer.)

Experimentally, with the following code

from sympy import *
var('x')
var('n', integer = True)

f = ( (2+(x+3) * sqrt(-x))/(2*(x+4)) *
      sqrt(-x)**n + ((2-(x+3) * sqrt(-x)) / (2*(x+4))) * (-sqrt(-x))**n +
      + (x+2+sqrt(x*(x+4)))/(2*(x+4)) *
        ( (x+sqrt(x*(x+4)))/2 )**n +
        (x+2-sqrt(x*(x+4)))/(2*(x+4)) *
        ( (x-sqrt(x*(x+4)))/2 )**n
        )
pprint(f)
for i in range(10):
   eq = simplify(f.subs(n,i))
   print ('========== n = ', i)
   pprint(eq)
   print ('have solutions')
   sols = solve(eq)
   pprint (sols)
   pprint ('approx. = ')
   pprint ([s.evalf() for s in sols])

we see the following: $$ \begin{array}{l} f_{0} = 1\\ f_{1} = 0\\ f_{2} = x^{2}\\ f_{3} = x^{2} \left(x + 2\right)\\ f_{4} = x^{2} \left(x^{2} + 2 x + 1\right)\\ f_{5} = x^{3} \left(x^{2} + 3 x + 1\right)\\ f_{6} = x^{4} \left(x^{2} + 4 x + 4\right)\\ f_{7} = x^{4} \left(x^{3} + 5 x^{2} + 7 x + 3\right)\\ f_{8} = x^{4} \left(x^{4} + 6 x^{3} + 11 x^{2} + 6 x + 1\right)\\ f_{9} = x^{5} \left(x^{4} + 7 x^{3} + 16 x^{2} + 13 x + 2\right)\\ \end{array} $$

And it seems that already $f_9$ have some imaginary roots, albeit very small. Do these imaginary roots really converge to 0 when $n\to \infty$ or they stay on the same magnitude ?