%I #28 Jan 06 2020 17:29:33

%S 1,6,84,1680,39240,999216,26899896,752939424,21691531800,638947312080,

%T 19155738105504,582589712312064,17930566188602136,557417298916695600,

%U 17477836958370383280,552090876791399769600,17552554240486710112920,561230779055361080132880

%N Binomial transform of the Apéry numbers (A005259).

%C Starting with the a(3) term, each term is divisible by 8. (Empirical observation.)

%C The above is true and follows easily from the pair of known congruences for the Apéry numbers A(n): A(2*n) == 1 (mod 8) and A(2n+1) == 5 (mod 8). - _Peter Bala_, Jan 06 2020

%H Jackson Earles, Justin Ford, Poramate Nakkirt, Marlo Terr, Dr. Ilia Mishev, Sarah Arpin, <a href="https://www.colorado.edu/math/binomial-transforms-sequences-fall-2018">Binomial Transforms of Sequences</a>, Fall 2018.

%H N. J. A. Sloane, <a href="/transforms.txt">Transforms</a>

%F a(n) ~ 2^(n - 3/4) * 3^(n + 3/2) * (1 + sqrt(2))^(2*n - 1) / (Pi*n)^(3/2). - _Vaclav Kotesovec_, Dec 17 2018

%F The Gauss congruences hold: a(n*p^k) == a(n*p^(k-1)) (mod p^k) for all primes p and n a positive integer. - _Peter Bala_, Jan 06 2020

%e a(2) = binomial(2,0)*A(0) + binomial(2,1)*A(1) + binomial(2,2)*A(2), where A(k) denotes the k-th Apéry number. Using this definition:

%e a(2) = binomial(2,0)*(binomial(0,0)*binomial(0,0))^2 + binomial(2,1)*((binomial(1,0)*binomial(1,0))^2 + (binomial(1,1)*binomial(2,1))^2) + binomial(2,2)*((binomial(2,0)*binomial(2,0))^2 + (binomial(2,1)*binomial(3,1))^2 + (binomial(2,2)*binomial(4,2))^2) = 84.

%t a[n_] := Sum[Binomial[n, k] * Sum[(Binomial[k, j] * Binomial[k+j, j])^2, {j, 0, k}], {k, 0, n}]; Array[a, 20, 0] (* _Amiram Eldar_, Dec 13 2018 *)

%o (Sage)

%o def OEISbinomial_transform(N, seq):

%o BT = [seq[0]]

%o k = 1

%o while k< N:

%o next = 0

%o j = 0

%o while j <=k:

%o next = next + ((binomial(k,j))*seq[j])

%o j = j+1

%o BT.append(next)

%o k = k+1

%o return BT

%o Apery = oeis('A005259')

%o OEISBinom = OEISbinomial_transform(18,Apery.first_terms(20))

%o (Julia)

%o function BinomialTransform(seq)

%o N = length(seq)

%o bt = Array{BigInt,1}(undef,N)

%o bt[1] = seq[1]

%o for k in 1:N-1

%o next = BigInt(0)

%o for j in 0:k next += binomial(k, j)*seq[j+1] end

%o bt[k+1] = next

%o end

%o bt end

%o BinomialTransform([A005259(n) for n in 0:18]) |> println # _Peter Luschny_, Jan 06 2020

%Y Cf. A005259, A322519.

%K nonn,easy

%O 0,2

%A _Sarah Arpin_, Dec 13 2018