inner_product#

skfda.misc.inner_product(arg1: Vector, arg2: Vector, *, _matrix: bool = False, _domain_range: Tuple[Tuple[float, float], ...] | None = None, **kwargs: Any) numpy.ndarray[Any, numpy.dtype[numpy.float64]]#

Return the usual (\(L_2\)) inner product.

Calculates the inner product between matching samples in two FDataGrid objects.

For two samples x and y the inner product is defined as:

\[\langle x, y \rangle = \sum_i x_i y_i\]

for multivariate data and

\[\langle x, y \rangle = \int_a^b x(t)y(t)dt\]

for functional data.

The two arguments must have the same number of samples, or one should contain only one sample (and will be broadcasted).

Parameters:
  • arg1 (Vector) – First sample.

  • arg2 (Vector) – Second sample.

  • kwargs (Any) – Additional parameters for the function.

  • _matrix (bool) –

  • _domain_range (Tuple[Tuple[float, float], ...] | None) –

Returns:

Vector with the inner products of each pair of samples.

Return type:

ndarray[Any, dtype[float64]]

Examples

This function can compute the multivariate inner product.

>>> import numpy as np
>>> from skfda.misc import inner_product
>>>
>>> array1 = np.array([1, 2, 3])
>>> array2 = np.array([4, 5, 6])
>>> inner_product(array1, array2)
32

If the arrays contain more than one sample

>>> array1 = np.array([[1, 2, 3], [2, 3, 4]])
>>> array2 = np.array([[4, 5, 6], [1, 1, 1]])
>>> inner_product(array1, array2)
array([32, 9])

The inner product of the \(f(x) = x\) and the constant \(y=1\) defined over the interval \([0,1]\) is the area of the triangle delimited by the the lines \(y = 0\), \(x = 1\) and \(y = x\), that is, \(0.5\).

>>> import skfda
>>>
>>> x = np.linspace(0,1,1000)
>>>
>>> fd1 = skfda.FDataGrid(x,x)
>>> fd2 = skfda.FDataGrid(np.ones(len(x)),x)
>>> inner_product(fd1, fd2)
array([ 0.5])

If the FDataGrid object contains more than one sample

>>> fd1 = skfda.FDataGrid([x, np.ones(len(x))], x)
>>> fd2 = skfda.FDataGrid([np.ones(len(x)), x] ,x)
>>> inner_product(fd1, fd2).round(2)
array([ 0.5, 0.5])

If one argument contains only one sample it is broadcasted.

>>> fd1 = skfda.FDataGrid([x, np.ones(len(x))], x)
>>> fd2 = skfda.FDataGrid([np.ones(len(x))] ,x)
>>> inner_product(fd1, fd2).round(2)
array([ 0.5, 1. ])

It also work with basis objects

>>> basis = skfda.representation.basis.MonomialBasis(n_basis=3)
>>>
>>> fd1 = skfda.FDataBasis(basis, [0, 1, 0])
>>> fd2 = skfda.FDataBasis(basis, [1, 0, 0])
>>> inner_product(fd1, fd2)
array([ 0.5])
>>> basis = skfda.representation.basis.MonomialBasis(n_basis=3)
>>>
>>> fd1 = skfda.FDataBasis(basis, [[0, 1, 0], [0, 0, 1]])
>>> fd2 = skfda.FDataBasis(basis, [1, 0, 0])
>>> inner_product(fd1, fd2)
array([ 0.5       , 0.33333333])
>>> basis = skfda.representation.basis.MonomialBasis(n_basis=3)
>>>
>>> fd1 = skfda.FDataBasis(basis, [[0, 1, 0], [0, 0, 1]])
>>> fd2 = skfda.FDataBasis(basis, [[1, 0, 0], [0, 1, 0]])
>>> inner_product(fd1, fd2)
array([ 0.5 , 0.25])