This site is supported by donations to The OEIS Foundation.

User:Peter Luschny/SequenceTransformations

From OeisWiki
Jump to: navigation, search

Transformations of
Integer Sequences

This article reconsiders and reimplements parts of

TRANSFORMATIONS OF INTEGER SEQUENCES
by N. J. A. Sloane (Maple version) and
Olivier Gerard (Mathematica version).

Note that some important transformations are defined differently here. This is always indicated in the text. The computer language used is the Sage dialect of Python.

Some sequence tools

Bisection

###########################################################
# name:   BISECTION
# param:  A sequence
# param:  k in {0,1}
# return: the k-th bisection of A
###########################################################

def bi_section(A, k) :
    L = []
    b = k == 0
    for a in A :
        if b : L.append(a)
        b = not b
    return L

Modulus-section

###########################################################
# name:   MODSECTION
# param:  A sequence
# param:  modulo m ≥ 2
# param:  k in {0,1,..,m-1}
# return: the k-th mod-section of A
###########################################################

def mod_section(A, m, k) :
    L = []
    b = 0
    for a in A :
        if b % m == k : L.append(a)
        b = b + 1
    return L

Coefficient list

##################################################################
# name:   COEFFICIENTLIST
# param:  p polynomial in the x variable with integer coefficients.
# return: the sequence of all the coefficients
##################################################################

def coefficient_list(p, x) :
    c = 0
    L = []
    R = p.coefficients(x)
    for r in R :
        while r[1] <> c :
            L.append(0)
            c = c + 1
        L.append(r[0])
        c = c + 1
    return L

Trans. based on arithmetic

First differences

###########################################################
# name:   FIRSTDIFFERENCE
# param:  A sequence
# return: The first differences of A
###########################################################

def first_difference(A) :
    if A == [] : return []
    L = []
    b = A[0]
    p = False
    for a in A :
        if p :
            L.append(a - b)
            b = a
        else : p = True
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 0,0,0,0,0,0,0,... 1,1,1,1,1,1,1,... A001563

Partial sums

###########################################################
# name:   PARTIALSUM
# param:  A sequence
# return: sequence of partial sums b[n] = a[0]+..+a[n].
###########################################################

def partial_sum(A) :
    if A == [] : return [0]
    p = 0
    L = []
    for a in A :
        p = a + p
        L.append(p)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 1,2,3,4,5,6,7,... A000217 A003422

Partial products

###########################################################
# name:   PARTIALPRODUCT
# param:  A sequence
# return: sequence of partial products b[n] = a[0]*..*a[n].
###########################################################

def partial_product(A) :
    if A == [] : return [1]
    p = 1
    L = []
    for a in A :
       p = p * a
       L.append(p)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 1,1,1,1,1,1,1,... 1,2,6,24,120,... A000178

partial_sum and partial_product are the prototypes of the following blueprint: Given an operation '#' and a special value X chosen by convention in the case of an 'empty operation'. Then 'partial operation #' applied to A = [a, b, c, ...] is expected to return:

if A = [] then [X] else partial#(A) = [a, a#b, a#b#c, ...].

In Sloane's transformation library PSUM([a, b, c]) returns [a, a+b, a+b+c] and behaves as expected whereas PRODS([a, b, c]) unexpectedly returns [1, a, a*b] (and thus is independent from the c in the input!).

Partial signed sums

###########################################################
# name:   PARTIALSUMSIGN
# param:  A sequence.
# return: sequence of
###########################################################

def partial_signsum(A) :
    if A == [] : return [0]
    p = 0
    L = []
    for a in A :
       p = a - p
       L.append(p)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000035 A008619 A058006

Partial alternating sums

###########################################################
# name:   ALTERNATINGSUM
# param:  A sequence.
# return: sequence of alternating partial sums,
#         b[n] = a[0]-a[1]+..+/-a[n].
###########################################################

def alternating_sum(A) :
    if A == [] : return [0]
    p = 0
    b = True
    L = []
    for a in A :
        if b : p = a + p
        else : p = p - a
        L.append(p)
        b = not b
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) ± A000035 ± A008619 ± A058006


Trans. based on connection coefficients

Binomial

###########################################################
# name:   BINOMIALTRANS
# param:  A sequence.
# return: the binomial transform of A
###########################################################

def binomial_trans(A) :
    L = []
    for n in range(len(A)) :
        s = sum(binomial(n,k)*A[k] for k in (0..n)) 
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000079 A001792 A000522

Binomial invers

###########################################################
# name:   BINOMIALINVTRANS
# param:  A sequence.
# return: the inverse binomial transform of A
###########################################################

def binomial_invtrans(A) :
    L = []
    for n in range(len(A)) :
        s = sum((-1)^(n-k)*binomial(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 1,0,0,0,0,0,0,... 1,1,0,0,0,0,... A000166

binomial_trans and binomial_invtrans are the prototypes of transformations based on integer triangles. Given A = a_{n, k} and B = b_{n, k} for n, k=0...L such that the matrices A and B are inverse to each other and a_{n, k} = b_{n, k} = 0 for all n < k ≤ L, then the associated transformation and it's inverse are described by

  u_n = sum a_{n,k} v_k  <==>  v_n = sum b_{n,k} u_k.

Based on this scheme two other transformations are easily given: the Stirling transformations and their inverses and the Lah transformation and it's inverse.

Stirling2

###########################################################
# name:   STIRLING2TRANS
# param:  A sequence.
# return: the Stirling2 transform of A
###########################################################

def stirling2_trans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum(stirling_number2(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000110 ~ A000110 A000670

Stirling2 invers

###########################################################
# name:   STIRLING2INVTRANS
# param:  A sequence.
# return: the inverse Stirling2 transform of A
###########################################################

def stirling2_invtrans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum((-1)^(n-k)*stirling_number1(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 1,1,0,0,0,0,0,... ~ A133942 A006252

stirling_number1 denotes here the unsigned version of these numbers, the Stirling cycles numbers.

Stirling1

###########################################################
# name:   STIRLING1TRANS
# param:  A sequence.
# return: the Stirling1 transform of A
###########################################################

def stirling1_trans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum(stirling_number1(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000142 A000774 A007840

Stirling1 invers

###########################################################
# name:   STIRLING1INVTRANS
# param:  A sequence.
# return: the inverse Stirling1 transform of A
###########################################################

def stirling1_invtrans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum((-1)^(n-k)*stirling_number2(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000587 1,2,1,-3,0,13,... A000012

Unfortunately the transformations STIRLING and STIRLINGi in Sloane's transformation library do not agree with the definitions above. Instead they compute the following transformations which have a shift in the indices of the Stirling numbers:

STIRLING

###########################################################
# name:   STIRLING
# param:  A sequence.
# return: Sloane's Stirling transform of A
###########################################################

def STIRLING(A) :
    L = []
    for n in range(len(A)) : 
        s = sum(stirling_number2(n+1,k+1)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000110 A005493 A000629

STIRLINGi

###########################################################
# name:   STIRLINGi
# param:  A sequence.
# return: Sloane's inverse Stirling transform of A
###########################################################

def STIRLINGi(A) :
    L = []
    for n in range(len(A)) : 
        s = sum((-1)^(n-k)*stirling_number1(n+1,k+1)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) 1,0,0,0,0,0,0,... ~ A000142 A089064


Lah

###########################################################
# name:   LAHTRANS
# param:  A sequence.
# return: the Lah transform of A
###########################################################
def lah_number(n, k) :  # defined for n ≥ 0, k ≥ 0.
    return (-1)^n*factorial(n-k)*binomial(n,n-k)*binomial(n-1,n-k)
1
0, -1
0, 2, 1
0, -6, -6, -1
0, 24, 36, 12, 1
def lah_trans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum(lah_number(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) ± A000262 A002720 ± A002866
T(T(A)) 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...

Lah invers

###########################################################
# name:   LAHINVTRANS
# param:  A sequence.
# return: the inverse Lah transform of A
###########################################################

def lah_invtrans(A) :
    L = []
    for n in range(len(A)) : 
        s = sum((-1)^(n-k)*lah_number(n,k)*A[k] for k in (0..n))
        L.append(s)
    return L
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) ± A066668 A202410 A154955
T(T(A)) 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...

By the way, the equation lah_invtrans([1,2,3,..]) = A202410 can be rewritten in mathematical terms, where Ln(x) denotes the n-th Laguerre polynomial, as:

This is an example how sequence transformations can be used to find new formulas.

Let's come back to the Lah transformations. What makes these transformations special is the involutory character of the Lah numbers. This is reflected by the transformations in the equation T(T(A)) = A.

For instance if we start with 1,2,3,.. and apply the Lah transformation we get the number of matchings in the complete bipartite graph K(n,n). If we apply the Lah transformation on this sequence again we are back at 1,2,3,.. Similarly if we apply the inverse Lah transformation to 1,2,3,.. we get n!(Ln(1)-2Ln-1(1)). Applying the inverse Lah transformation on this sequence we are back at 1,2,3,.. .

Unfortunately Sloane chose in his transformation library to call different transformations the Lah resp. inverse Lah transformation. This naming must be rejected as these transformations do not possess the involution property characterizing Lah transformations.

Trans. based on generating functions

Euler

###########################################################
# name:   EULERTRANS
# param:  A sequence.
# return: the Euler transform of A
###########################################################

def euler_trans(A) :
    L = []; M = []
    for i in range(len(A)) : 
        s = sum(d*A[d-1] for d in divisors(i+1))
        L.append(s)
        s = s + sum(L[j-1]*M[i-j] for j in (1..i)) 
        M.append(s/(i+1))
    return M
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000041 A000219 A179327

Euler invers

###########################################################
# name:   EULERINVTRANS
# param:  A sequence.
# return: the inverse Euler transform of A
###########################################################

def euler_invtrans(A) :
    L = []; M = []
    for i in range(len(A)) :
        s = (i+1)*A[i] - sum(L[j-1]*A[i-j] for j in (1..i)) 
        L.append(s)
        s = sum(moebius((i+1)/d)*L[d-1] for d in divisors(i+1))
        M.append(s/(i+1))
    return M
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000007 1,1,1,0,0,-1,0,0,0,... 1,0,1,4,18,95,...

A nice application of the Euler and inverse Euler transformation is Gordon's theorem, a generalization of the Rogers-Ramanujan identities. Observe that euler_invtrans(A035937) generates the list [0,1,1,1,1,0,0], in endless repetition. Similarly euler_invtrans(A035939) generates the list [1,1,0,0,1,1,0] in endless repetition. This is just a way Gordon's theorem presents itself in some special cases. Conversely euler_trans can be used to generate the sequences A035937 through A035941 by applying some short 0-1 lists repeatedly.

def gordons_theorem(A, n) :
    L = []; M = [];
    m = len(A)
    for i in range(n) :
        s = sum(d*A[(d-1) % m] for d in divisors(i+1))
        L.append(s)
        s = s + sum(L[d-1]*M[i-d] for d in (1..i))
        M.append(s/(i+1))
    return M
    
def A035937_list(len) : return gordons_theorem([0,1,1,1,1,0,0], len)
def A035938_list(len) : return gordons_theorem([1,0,1,1,0,1,0], len)
def A035939_list(len) : return gordons_theorem([1,1,0,0,1,1,0], len)
def A035940_list(len) : return gordons_theorem([0,1,1,1,1,1,1,0,0], len)
def A035941_list(len) : return gordons_theorem([1,0,1,1,1,1,0,1,0], len)

Comparing euler_trans with gordons_theorem we see that the only difference is the build-in repetition modulo the length of the generating 0-1 list.

        s = sum(d*A[d-1] for d in divisors(i+1))   
        s = sum(d*A[(d-1) % m] for d in divisors(i+1))

Finally let's look at three simple 0-1 'generators' in the above sense: euler_trans([1,0,1,0,1,0,1,..]) = number of partitions of n into odd parts, A000009, and euler_trans([0,1,0,1,0,1,0,..]) = number of partitions of n into even parts, A035363. So it is no surprise now that euler_trans([1,1,1,1,1,...]) is the number of partitions of n, A000041. But it is a surprise that these simple characterizations are not mentioned in the database.

The transformations motzkin_trans, catalan_trans and their inverses are based on Richard J. Mathar's Maple programs 'transforms3'.

Motzkin

###########################################################
# name:   MOTZKINTRANS
# param:  A sequence.
# return: the Motzkin transform of A
###########################################################
def motzkin_trans(A) :
    d = 0
    p = 0  
    l = len(A) - 1
    f = (1 - x - sqrt(1 - 2*x - 3*x^2))/2/x
    m = f.taylor(x, 0, l)
    
    for a in A : 
        p += (a*m^d).taylor(x, 0, l)
        d += 1
    return coefficient_list(p, x)
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A)  A005773  A036908 1,1,3,12,56,303, ...

Motzkin invers

###########################################################
# name:   MOTZKININVTRANS
# param:  A sequence.
# return: the inverse Motzkin transform of A
###########################################################
def motzkin_invtrans(A) :
    d = 0
    p = 0  
    l = len(A) - 1
    f = x / (1 + x + x^2)
    m = f.taylor(x, 0, l)
    
    for a in A : 
        p += (a*m^d).taylor(x, 0, l)
        d += 1
    return coefficient_list(p, x)
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A117569 A147657 1,1,1,2,9,45,268, ..

Catalan

###########################################################
# name:   CATALANTRANS
# param:  A sequence.
# return: the Catalan transform of A
###########################################################
def catalan_trans(A) :
    d = 0
    p = 0  
    l = len(A) - 1
    f = (1 - sqrt(1 - 4*x))/2
    m = f.taylor(x, 0, l)
    
    for a in A : 
        p += (a*m^d).taylor(x, 0, l)
        d += 1
    return coefficient_list(p, x)
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A000108 A000108  1,1,3,12,57,312,...

Catalan inverse

###########################################################
# name:   CATALANINVTRANS
# param:  A sequence.
# return: the inverse Catalan transform of A
###########################################################
def catalan_invtrans(A) :
    d = 0
    p = 0  
    l = len(A) - 1
    f = x*(1 - x)
    m = f.taylor(x, 0, l)
    
    for a in A : 
        p += (a*m^d).taylor(x, 0, l)
        d += 1
    return coefficient_list(p, x)
A 1,1,1,1,1,1,1,... 1,2,3,4,5,6,7,... 1,1,2,6,24,120,...
T(A) A010892 A099254 A013999

Sequence operations

Convolution

###########################################################
# name:   CONVOLUTION
# param:  A, B sequences.
# return: the convolution of A and B
###########################################################
def convolution(A, B) :
    m = min(len(A), len(B))
    C = []
    for n in range(m) :
        s = sum(A[k]*B[n-k] for k in (0..n))
        C.append(s)
    return C
A  0,1,2,3,4,5,6,7,...  A002802 A000035
B  1,1,2,6,24,120,... A000984 1,2,3,4,5,6,7,...
P(A, B) A014144   A038845   A035048  

Binomial convolution

###########################################################
# name:   BINOMIALCONVOLUTION
# param:  A, B sequences.
# return: the convolution of A and B
###########################################################
def binomial_convolution(A, B) :
    m = min(len(A), len(B))
    C = []
    for n in range(m) :
        s = sum(binomial(n,k)*A[k]*B[n-k] for k in (0..n))
        C.append(s)
    return C

If B[n] = 1 for all n, then the binomial convolution reduces to the binomial transformation.

A  0,1,2,3,4,5,6,...  1,1,1,1,1,1,1,... 1,2,3,4,5,6,...
B  1,1,2,6,24,...  1,1,2,6,24,...  1,1,2,6,24...
P(A, B) A007526   A000522  A111063  


Divconv numerator

###########################################################
# name:   DIVCONVNUMERATOR
# param:  A, B sequences, B[n] <> 0 for all n.
# return: numerator of sum(A[k]/B[n-k],k=0..n)
###########################################################
def divconv_numerator(A, B) :
    m = min(len(A), len(B))
    C = []
    for n in range(m) :
        s = sum(A[k]/B[n-k] for k in (0..n))
        C.append(s.numerator())
    return C
A  1,2,3,4,5,6,...  1,1,1,1,1,1,1,...  1,1,1,1,1,1,1,...
B  1,2,3,4,5,6,...  1,2,3,4,5,6,7,...  1,1,2,6,24,...
P(A, B) A027612 A001008  A061354

Divconv denominator

###########################################################
# name:   DIVCONVDENOMINATOR
# param:  A, B sequences, B[n] <> 0 for all n.
# return: denominator of sum(A[k]/B[n-k],k=0..n)
###########################################################
def divconv_denominator(A, B) :
    m = min(len(A), len(B))
    C = []
    for n in range(m) :
        s = sum(A[k]/B[n-k] for k in (0..n))
        C.append(s.denominator())
    return C
A  1,2,3,4,5,6,...  1,1,1,1,1,1,1,...  1,1,1,1,1,1,1,...
B  1,2,3,4,5,6,...  1,1,2,6,24,...  1,1,2,6,24,...
P(A, B) A027611 A002805  A061355

Sequence to triangle transformations

Up to now all transformations had the form
     T ← trans(A),  where A is a list A = [a0, a1, ..., am] or
     T ← trans(A, B),  where B is the list B = [b0, b1, ..., bm].

However, this format is quite restrictive. Our next example uses a more general setup: a sequence of polynomials pn with degree(pn) = n which have coefficients  depending on a given sequence f (which we will implement this time as a function).

This gives much more flexibility; in particular one can look now at

  •  the polynomials themselves,
  •  the coefficients of the polynomials (which will give triangles) and 
  •  the list of values of the polynomials evaluated at some specific argument.

Our first example are the A-polynomials (this is just an ad-hoc name). Corresponding to the three views on the transformation we implement three functions.

A-polynomials

###########################################################
# name:   APOLYNOMIALS  (ad hoc name)
# param:  f function, n integer, x variable.
###########################################################
def A_polynomial(f, n, x) :
    s = add(add((-1)^v*binomial(k, v)*f(k+1)*(x+v+1)^n 
    for v in (0..k)) for k in (0..n))
    return expand(s) 
    
def A_triangle(f, n) :
    var('x')
    return coefficient_list(A_polynomial(f, n, x), x)

def A_list(f, len, a) :
    return [A_polynomial(f, n, x=a) for n in (0..len)]
f1
f1X+f1−f2
f1X2+(2f1−2f2)X+2f3+f1−3f2
f1X3+(3f1−3f2)X2+(6f3+3f1−9f2)X−6f4+12f3+f1−7f2
f1X4+(4f1−4f2)X3+(12f3+6f1−18f2)X2 +(−24f4+48f3+4f1−28f2)X−60f4+50f3+f1−15f2+24f5

Some examples:

  • Let f generate the constant sequence 1,1,1,1,...
    def uno(n) : return 1.
    • A_polynomial(uno, n, x) is 1, x, x^2, x^3,...
    • A_triangle(uno, n) is (as a matrix) the identity matrix.
    • A_list(uno, len, n) is the sequence n^k, k=0,1,2,..
  • Let f generate the sequence 1,2,3,4,...
    def id(n) : return n.
    • A_triangle(id, n) is A130595
    • A_list(id, len, -n) is the sequence (-n-1)^k, k=0,1,2,..
  • Let f generate the harmonic numbers,
    def H(n) : return add(1/i for i in (1..n)).
    • A_polynomial(H, n, x) Lo and behold! These are the Bernoulli polynomials!
    • The numerators of A_list(H,len,-1) is A164558 with different signs.
  • Let f = skp, defined as
    def skp(m) : if m % 4 == 0 : return 0 else : return 1/((-1)^(m//4)*2^((m-1)//2))
    • A_polynomial(skp, n, x) are the Swiss knife polynomials.
    • A_triangle(skp, n) is A153641 (zeros to be filled in).
    • A_list(skp, len, 0) are the Euler numbers .
    • A_list(skp, len, 1) are the tangent numbers A155585, A000182.

The next transformation should also be in the quiver of every IS-researcher. I do not know a name or reference at present for this transformation. I call it here c_trans but a better name certainly is desirable. Please let me know if you know a reference.

C-transformation

###########################################################
# name:   CTRANSFORMATION  (another ad hoc name)
# param:  A sequence
# return: the transform of A which is a triangle
###########################################################
def c_trans(A) :    
    L = len(A) + 1
    C = [1 for i in (0..L)]; C[0] = 0
    
    for k in (1..L) :
        for n in range(k-1,0,-1) :
            C[n] = C[n-1] + C[n+1]*A[n-1]
 
        print [C[n] for n in (1..k)]
A 1,1,1,1,1,... 1,2,3,4,5,6,... A000035 A059841
T(A) A033184    A001497    A103631     A065941  

Deleham delta

###########################################################
# name:   DELEHAMDELTA
# param:  r, s sequences
# return: r Δ s, a triangle
###########################################################
def deleham_delta(R, S) :
    L = min(len(R),len(S)) + 1
    A = [SR(R[k] + x*S[k]) for k in range(L-1)]
    C = [SR(1) for i in range(L+1)]; C[0] = SR(0)
    
    for k in (1..L) :
        for n in range(k-1,0,-1) :
            C[n] = C[n-1] + C[n+1]*A[n-1]
 
        p = expand(C[1])
        print [p.coefficient(x,n) for n in (0..k-1)]
A  1,1,1,1,1,1,1,... A110654 A000035
B  1,1,1,1,1,1,1,... A000007 A059841
A Δ B A085880 A084938 A090181

Deleham defines an operation mapping two sequences to a triangle. He first encodes the two given sequences R and S into a polynomial A,
  A[k] = R[k] + x*S[k],
then applies c_trans to A (and gets by that a triangle of polynomials).  Finally he decodes the polynomials in the first column by extracting their coefficients.

(For readers who want to read the above code but are not fluent in Sage parlance: The expression SR(R[k] + x*S[k]) means: 'Interpret R[k]+x*S[k] as a polynomial (in the Symbolic Ring)'. In Maple for example such a coercion is not necessary, one would just write R[k]+x*S[k].)

Two remarks are in order. First, our description is somewhat simpler than Deleham's original description, which he gave in A084938. A two dimensional array is not necessary to implement the recursion; also the polynomial here is not in two variables x and y but only in x, simplifying things further.

Second, the c_trans applied to R[k] + x*S[k] generates a triangle of polynomials, thus the natural way to transform it back to integers would be to extract the coefficients of all these polynomials. Deleham selects only one column in the triangle, which might be right choice in many circumstances; however it is not imperative to do so and other choices might also be interesting.

With Maple:

DELEHAMDELTA := proc(r, s) local L, k, A, p, R, C, n, W;
R := NULL; C[0] := 0; L := min(nops(r), nops(s))+1;

for k from 0 to L-2 do
    A[k] := x*r[k+1] + s[k+1];
od:

for k from 1 to L do
    C[k] := 1;
    for n from k-1 by -1 to 1 do
        C[n] := C[n-1] + C[n+1]*A[n-1]
    od;

    p := expand(C[1]);
    W := seq(coeff(p,x,k-n), n=1..k);
    R := R, W;
    lprint(W);
od:
R end:

Miscellaneous transformations

Exponential transform

###########################################################
# name:   EXP-TRANSFORM
# param:  s sequence
# return: the transform of s 
###########################################################

The exponential transform is defined as (Maple)

exptrans := proc(s) local g; g := proc(n) option remember; local k;
`if`(n=0, 1, add(binomial(n-1, k-1)*s(k)*g(n-k), k=1..n)) end end:

Some sequences which can be generated with this transformation are: A000110, the Bell numbers (number of set partitions of {1..n}); A000587, the Rényi numbers (number of set partitions of {1..n} with an even number of parts, minus the number of partitions with an odd number of parts); A000248, the number of forests with n nodes and height at most 1; A007837, the number of partitions of n-set with distinct block sizes; and A005651, the sum of multinomial coefficients. 

A000110 := n-> exptrans(m->1)(n):  
A000587 := n-> exptrans(m->-1)(n): 
A000248 := n-> exptrans(m->m)(n): 
A007837 := n-> exptrans(m->A182927(m))(n): 
A005651 := n-> exptrans(m->A182926(m))(n):
A 1,1,1,1,1,... -1,-1,-1,-1,... 0,1,2,3,4,... A182927 A182926
T(A) A000110    A000587    A000248     A007837    A005651  

If the reader is in an explorative mood he might also wish to find a combinatorial interpretation for

A_new := n-> exptrans(m->A038041(m))(n): 
1, 1, 3, 9, 38, 168, 915, 5225, 34228, 236622, 1805297,...

Iterated exponential transform

###########################################################
# name:   ITERATED EXP-TRANSFORM
# param:  s sequence
# return: the transform of s 
###########################################################

The iterated exponential transform is the exponential transform applied n times to s and then evaluated at k. 

iterexp := s -> proc(n,k) (exptr@@n)(m->s(m))(k) end:

Sequences which can be generated with this transformation are the iterated Bell numbers A144150 and the iterated Rényi numbers A210638.

A144150 := iterexp(m-> 1): seq(lprint(seq(A144150(k,n),k=0..6)),n=0..6);
A210638 := iterexp(m->-1): seq(lprint(seq(A210638(n,k),k=0..6)),n=0..6);
A209631 := iterexp(m-> m): seq(lprint(seq(A209631(n,k),k=0..6)),n=0..6);
A 1,1,1,1,1,... -1,-1,-1,-1,... 0,1,2,3,4,...
T(A) A144150    A210638    A209631   


The Hankel transformation which we will consider here was not first introduced by Hermann Hankel (Hankel introduced a different transform) but rather "by the author [John W. Layman] in sequence A055878 of the On-Line Encyclopedia of Integer Sequences" according to this article. Here we introduce three variants of the transformation of integer sequences.

Hankel

###########################################################
# name:   HANKEL
# param:  f a sequence generator
# n:      length of transform
# return: Hankel integer transform of [f(0)..f(n-1)]
###########################################################

def hankel_even(f, n) :
    R = []
    for k in range(n) :
        R.append(matrix(k, k, lambda i, j: f(i+j)).det())
    return R
    
def hankel_odd(f, n) :
    R = []
    for k in range(n) :
        R.append(matrix(k, k, lambda i, j: f(i+j+1)).det())
    return R    
    
def hankel(f, n) :
    R = []
    for k in range(n) :
        R.append(matrix(k, k, lambda i, j: f(i+j)).det())
        R.append(matrix(k, k, lambda i, j: f(i+j+1)).det())
    return R

It is clear that hankel_even and hankel_odd are just the bisections of hankel. For example 

def catalan(n) : return binomial(2*n, n)/(n+1)

print hankel_even(catalan,9)
print hankel_even(factorial,6)
    [1, 1, 1, 1, 1, 1, 1, 1, 1]
    [1, 1, 1, 4, 144, 82944]

print hankel_odd(catalan,9)
print hankel_odd(factorial,6)
    [1, 1, 1, 1, 1, 1, 1, 1, 1]
    [1, 1, 2, 24, 3456, 9953280]

print hankel(catalan,9)  
print hankel(factorial,6)
    [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
    [1, 1, 1, 1, 1, 2, 4, 24, 144, 3456, 82944, 9953280]

What Layman calls the Hankel transform in his paper and Sloane HANKEL in his Maple library is our hankel_even. The reason why we do not follow Layman/Sloane are the following facts:

  • The Layman/Sloane Hankel transformation is many-to-one.
  • The Hankel transformation as defined here uniquely determines the original number sequence provided that all determinants are not zero.

Thus our definition is much more useful in the study sequences. For instance the sequence of Catalan numbers is the unique sequence which maps to 1,1,1,1,... This is a beautiful and important characterization of the Catalan numbers.

The definition of Layman jumps too short; it captures only one side of the coin. The definition used here was (to the best of my knowledge) first introduced by Martin Aigner who generalized the Hankel integer transformation into a much larger framework.

Three more examples taken from Paul Barry's A120580 might illustrate what has been said.

def A006134(n) : 
    return sum(binomial(2*k, k) for k in (0..n))
def A098479(n) : 
    return sum(binomial(n-k,k)*binomial(n-2*k,k) for k in (0..(n//2)))
def A025565(n) : 
    return sum(binomial(n-1,k)*binomial(n-k+1,k+1) for k in (0..(n//2)))
    
print hankel_even(A006134,9)
print hankel_even(A098479,9)
print hankel_even(A025565,9)
    [1, 1, 0, -4, -8, 0, 32, 64, 0]
    [1, 1, 0, -4, -8, 0, 32, 64, 0]
    [1, 1, 0, -4, -8, 0, 32, 64, 0]

print hankel_odd(A006134,9)
print hankel_odd(A098479,9)
print hankel_odd(A025565,9)
    [1, 3, 6, -8, -72, -144, 64, 960, 1920]
    [1, 1, 2, -8, -24,  -48, 64, 320,  640]
    [1, 2, 4,  8,  16,   32, 64, 128,  256]

print hankel(A006134,9)
print hankel(A098479,9)
print hankel(A025565,9)
    [1,1,1,3,0,6,-4,-8,-8,-72,0,-144,32,64,64,960,0,1920]
    [1,1,1,1,0,2,-4,-8,-8,-24,0, -48,32,64,64,320,0, 640]
    [1,1,1,2,0,4,-4, 8,-8, 16,0,  32,32,64,64,128,0, 256]

Next month we will look into number-theoretic transformations.

Retrieved from "https://oeis.org/w/index.php?title=User:Peter_Luschny/SequenceTransformations&oldid=1570513"