Proving integrality of a sequence of numbers How would one go about proving whether or not the terms of a sequence such as the following
$$1, 1, 1, 1, 1, 6, 21, 181, 5221, 1090981, 986241401, 51490676208426,\ldots$$
defined by
$$a_na_{n-5}=a_{n-1}(a_{n-2}+a_{n-3}+a_{n-4})+a_{n-2}(a_{n-3}+a_{n-4})+a_{n-3}a_{n-4},$$
with $a_1=a_2=\ldots=a_5=1$, are always integers? I have checked the first 50 terms using Mathematica, but the numbers get huge pretty fast and my poor old pc cannot compute much more than that.
By the way, in the general case
$$a_na_{n-k}=a_{n-1}(a_{n-2}+\ldots+a_{n-k+1})+a_{n-2}(a_{n-3}+\ldots+a_{n-k+1})+\ldots+a_{n-k+2}a_{n-k+1},$$
with $a_1=a_2=\ldots=a_k=1$, it seems to be the case that sequences generated with even $k$ quickly get non-integer terms.
Realated: Somos sequences.
 A: From Section 7.6 of Alman, Cuenca and Huang, I extract the following proof without insight:
Define
$$b_m = \frac{a_m + a_{m+1} + a_{m+2} + a_{m+3} + a_{m+4}}{a_{m+1} a_{m+3}}.$$
A computer algebra system checks that $b_m = b_{m+2}$ (this is just an equality between two rational functions of $(a_m, a_{m+1}, a_{m+2}, a_{m+3}, a_{m+4})$.) In our case, $b_i$ alternates $5$, $10$, $5$, $10$, etcetera.
In particular, we have the recursion
$$a_{m+4} = (\mbox{$5$ or $10$}) \cdot a_{m+1} a_{m+3} - a_m - a_{m+1} - a_{m+2} - a_{m+3}$$
whose values are plainly integral.
A: This is Theorem 3.9.(6) of Alman, Cuenca and Huang with $n=5$ and $A=B=0$. I'm afraid I do not understand their proof well enough to say more.

I can say a little more to fill in context. Everything that I am saying is also in the introduction to Alman, Cuenca and Huang. A stronger property than integrality is Laurentness: Define rational functions $x_i$ by making $x_1$, $x_2$, $x_3$, $x_4$, $x_5$ formal variables and then setting
$$x_n = \frac{\sum_{1 \leq i < j \leq 4} x_{n-i} x_{n-j}}{x_{n-5}}.$$
Laurentness is the claim that all of these rational functions have denominators that are monomials in $(x_1, \ldots, x_5)$. In particular, we plug in $(x_1, x_2, x_3, x_4, x_5) = (1,1,1,1,1)$, then all the terms are integral
Fomin and Zelevinsky, in their "caterpillar lemma" describe a general recursive set up with the following properties: One starts with $((x_1, \ldots, x_n), (P_1, \ldots, P_n)$ where $P_i$ are polynomials in the $x_i$ obeying certain compatabilities. One then chooses an index $i$ to "mutate" at. One replaces $(x_1, \ldots, x_n)$ by $(x_1, \ldots, x_{i-1}, x'_i, x_{i+1}, \ldots, x_n)$ where
$$x'_i = \frac{P_i(x_1, \ldots, x_n)}{x_i}.$$
One also replaces $(P_1, \ldots, P_n)$ by modified polynomials according to a recipe which is too complex to give here. The caterpillar lemma then says that, no matter how many times you perform this operation, you will get $x$'s which are Laurent in the original $(x_1, \ldots, x_n)$. Lam and Pylyavskyy made certain technical improvements on Fomin and Zelevinsky, and ACH actually use the Lam-Pylyavskyy formalism. 
In particular, suppose that the FZ recipe, when mutating at $x_1$, turns $(P_1, P_2, \ldots, P_n)$ into $(P_n, P_1, P_2, \ldots, P_{n-1})$. One could then periodically mutate at $1$, $2$, $3$, ..., $n$, $1$, $2$, ..., $n$, ... You will then get a sequence $x_1$, $x_2$, ..., which obeys the Laurentness property, defined by the recursion
$$x_m = \frac{P_1(x_{m-n+1},\ldots, x_{m-1})}{x_{m-n}} \quad (\ast).$$
Conversely, given a Laurent sequence of the form $(\ast)$, one can wonder whether $P_1$ can be extended to an $n$-tuple $(P_1, \ldots, P_n)$ so the the FZ operation turns it into $(P_n, P_1, P_2, \ldots, P_{n-1})$. 
When this can be done, Alman, Cuenca and Huang say that $P_1$ is $1$-periodic. They give an algorithm to test whether a polynomial is $1$-periodic, and Theorem 3.9 lists several cases they have discovered using their algorithm.
