Non-negative low-rank approximations for
multi-dimensional arrays on statistical manifold
Kazu Ghalamkari1,2,

Mahito Sugiyama1,2

1 : The Graduate University for Advanced Studies, SOKENDAI
2 : National Institute of Informatics

International Conference on Information Geometry for Data Science (IG4DS 2022)

Motivation□ Non-negative low-rank approximation of data with various structures
Approximates with a linear combination of fewer bases (principal components) for feature
extraction, memory reduction, and pattern discovery.😀

≃

≃

≃

≃

2

Motivation□ Non-negative low-rank approximation of data with various structures
Approximates with a linear combination of fewer bases (principal components) for feature
extraction, memory reduction, and pattern discovery.😀

≃

≃

≃

≃
Non-negative constraint improves interpretability
3Motivation□ Non-negative low-rank approximation of data with various structures
Approximates with a linear combination of fewer bases (principal components) for feature
extraction, memory reduction, and pattern discovery.😀

≃

≃

≃

≃

Non-negative constraint improves interpretability
Low-rank approximation with non-negative constraints are based on gradient methods.
→ Appropriate settings for stopping criteria, learning rate, and initial values are necessary4 😢

Strategy□ Modeling with probability mass function on Directed Acyclic Graph(DAG).
□ Utilize projection theory of information geometry.

5Strategy□ Modeling with probability mass function on Directed Acyclic Graph(DAG).
□ Utilize projection theory of information geometry.

6

Strategy□ Modeling with probability mass function on Directed Acyclic Graph(DAG).
□ Utilize projection theory of information geometry.

7

Strategy□ Modeling with probability mass function on Directed Acyclic Graph(DAG).
□ Utilize projection theory of information geometry.

8

Contribution
Information Geometric Analysis using Distributions on DAGs that Correspond to Data Structures□ LTR: Faster Tucker-rank Reduction

No worries about initial values, stopping criterion and learning rate 😄

9Contribution
Information Geometric Analysis using Distributions on DAGs that Correspond to Data Structures□ LTR: Faster Tucker-rank Reduction

Rank-1 = rank 1,1,1

No worries about initial values, stopping criterion and learning rate 😄

10Contribution
Information Geometric Analysis using Distributions on DAGs that Correspond to Data Structures□ LTR: Faster Tucker-rank Reduction □ A1GM: Faster rank-1 missing NMF

Missing value

Find the most dominant factor rapidly.
Solve the task as a coupled NMF.
Rank-1 = rank 1,1,1

No worries about initial values, stopping criterion and learning rate 😄

11Contents□ Motivation, Strategy, and Contributions□ Introduction of log-linear model on DAG□ The best rank-1 approximation formula□ The best rank-1 NMMF

□ Legendre Tucker-Rank Reduction(LTR)□A1GM: faster rank-1 missing NMF

github.com/gkazunii/Legendre-tucker-rank-reduction

github.com/gkazunii/A1GM□ Theoretical Remarks□ Conclusion
12

3:00

Modeling tensor and matrix□ Flexible modeling is required to capture the structure of various data

Formulate low-rank approximations with probabilistic models on DAGs

13

Log-linear model on Directed Acyclic Graph (DAG)□DAG(poset)

is a DAG ⇔ for all 𝑠1 , 𝑠2 , 𝑠3 ∈the following three properties are satisfied.

(1) Reflexivity ∶ 𝑠1 ≤ 𝑠1 (2) Antisymmetry: 𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠1 ⇒ 𝑠1 = 𝑠2 (3)Transitivity:𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠3 ⇒ 𝑠1 ≤ 𝑠3

Mahito Sugiyama, Hiroyuki Nakahara and Koji Tsuda
14"Tensor balancing on statistical manifold“(2017)
ICML.

Log-linear model on Directed Acyclic Graph (DAG)□DAG(poset)

is a DAG ⇔ for all 𝑠1 , 𝑠2 , 𝑠3 ∈the following three properties are satisfied.

(1) Reflexivity ∶ 𝑠1 ≤ 𝑠1 (2) Antisymmetry: 𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠1 ⇒ 𝑠1 = 𝑠2 (3)Transitivity:𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠3 ⇒ 𝑠1 ≤ 𝑠3

□ log-linear model on DAG
We define the log-linear model on a DAG

as a mapping 𝑝:→ 0,1 ．Natural parameters 𝜽 describe the model.

𝜃-space

Mahito Sugiyama, Hiroyuki Nakahara and Koji Tsuda
15"Tensor balancing on statistical manifold“(2017)
ICML.

Log-linear model on Directed Acyclic Graph (DAG)□DAG(poset)

is a DAG ⇔ for all 𝑠1 , 𝑠2 , 𝑠3 ∈the following three properties are satisfied.

(1) Reflexivity ∶ 𝑠1 ≤ 𝑠1 (2) Antisymmetry: 𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠1 ⇒ 𝑠1 = 𝑠2 (3)Transitivity:𝑠1 ≤ 𝑠2 , 𝑠2 ≤ 𝑠3 ⇒ 𝑠1 ≤ 𝑠3

□ log-linear model on DAG
We define the log-linear model on a DAG

as a mapping 𝑝:→ 0,1 ．Natural parameters 𝜽 describe the model.

We can also describe the model by expectation parameters 𝜼 with Möbius function.

𝜃-space𝜂-space

Mahito Sugiyama, Hiroyuki Nakahara and Koji Tsuda
16"Tensor balancing on statistical manifold“(2017)
ICML.

Contents□ Motivation, Strategy, and Contributions□ Introduction of log-linear model on DAG□ The best rank-1 approximation formula□ The best rank-1 NMMF

□ Legendre Tucker-Rank Reduction(LTR)□A1GM: faster rank-1 missing NMF

github.com/gkazunii/Legendre-tucker-rank-reduction

github.com/gkazunii/A1GM□ Theoretical Remarks□ Conclusion
17

4:30

Introducing DAGs for Tensor

18

Introducing DAGs for Tensor

19Introducing DAGs for Tensor

20

Describe a tensor with (θ,η)

21

Describe a tensor with (θ,η)

22

Describe a tensor with (θ,η)

Möbius inversion formula23Describe a tensor with (θ,η)

Möbius inversion formula24Describe a tensor with (θ,η)

Relation between distribution and tensor

Möbius inversion formula

Random variables ： 𝑖, 𝑗, 𝑘 , indices of the tensor
Sample space： index set
Probability values ： tensor values 𝒫𝑖𝑗𝑘
25

One-body and many-body parameters

One-body parameter

Many-body parameter26

𝜽-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 subspace

One-body parameter

Many-body parameter27

𝜽-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 subspace
is e-flat. The projection is unique.
One-body parameter

Many-body parameter28

𝜽-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 subspace
is e-flat. The projection is unique.
One-body parameter

Many-body parameter

We can find the projection destination by a gradient-method.
But gradient-methods require Appropriate settings for stopping criteria, learning rate, and initial values 😢
29𝜽-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 subspace
is e-flat. The projection is unique.
One-body parameter

Many-body parameter

We can find the projection destination by a gradient-method.
But gradient-methods require Appropriate settings for stopping criteria, learning rate, and initial values 😢

Let us describe the rank-1 condition with the 𝜂-parameter.

30𝜼-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 condition (𝜼- representation)
Rank-1 subspace

𝜂𝑖𝑗𝑘 = 𝜂𝑖11 𝜂1𝑗1 𝜂11𝑘

One-body parameter

Many-body parameter31𝜼-representation of rank-1 tensor
Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

=

Rank-1 condition (𝜼- representation)=

𝜂𝑖𝑗𝑘 = 𝜂𝑖11 𝜂1𝑗1 𝜂11𝑘=

One-body parameter

Rank-1 subspace

Many-body parameter

The m-projection does not change one-body η-parameterShun-ichi Amari, Information Geometry and Its Applications, 2008, Theorem 11.6

32

Find the best rank-1 approximation

Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 condition (𝜼- representation)
Rank-1 subspace𝜂ҧ𝑖𝑗𝑘 = 𝜂ҧ𝑖11 𝜂ҧ1𝑗1 𝜂ҧ11𝑘

One-body parameter

Many-body parameter

33

Find the best rank-1 approximation

Rank-1 condition (𝜽-representation)
Its all many-body 𝜃-parameters are 0.

Rank-1 condition (𝜼- representation)
Rank-1 subspace𝜂ҧ𝑖𝑗𝑘 = 𝜂ҧ𝑖11 𝜂ҧ1𝑗1 𝜂ҧ11𝑘= 𝜂𝑖11 𝜂1𝑗1 𝜂11𝑘
One-body parameter

Many-body parameter
Möbius inversion formulaAll 𝜼-parameters after the projection are identified.
Using inversion formula, we found the projection destination.

34

Mean-field approximation and rank-1 approximation
Best rank-1 tensor formula for minimizing KL divergence (𝑑 = 3 )The best rank-1 approximation of 𝒫 ∈ ℝ𝐼×𝐽×𝐾
is given as>0which minimizes KL divergence from 𝒫.
We reproduce the result in K.Huang, et al. "Kullback-Leibler principal component for tensors is not NP-hard." ACSSC 2017

9:00
35

Mean-field approximation and rank-1 approximation
Best rank-1 tensor formula for minimizing KL divergence (𝑑 = 3 )The best rank-1 approximation of 𝒫 ∈ ℝ𝐼×𝐽×𝐾
is given as>0

By the way,
Frobenius error
minimization
is NP-hardwhich minimizes KL divergence from 𝒫.
We reproduce the result in K.Huang, et al. "Kullback-Leibler principal component for tensors is not NP-hard." ACSSC 201736

Mean-field approximation and rank-1 approximation
Best rank-1 tensor formula for minimizing KL divergence (𝑑 = 3 )The best rank-1 approximation of 𝒫 ∈ ℝ𝐼×𝐽×𝐾
is given as>0Normalized vector
depending on only 𝑖

Normalized vector
depending on only 𝑗

By the way,
Frobenius error
minimization
is NP-hard

Normalized vector
depending on only 𝑘which minimizes KL divergence from 𝒫.
We reproduce the result in K.Huang, et al. "Kullback-Leibler principal component for tensors is not NP-hard." ACSSC 2017

A tensor with 𝑑 indices is a joint distribution with 𝑑 random variables.
A vector with only 1 index is an independent distribution with only one random variable.

37

Mean-field approximation and rank-1 approximation
Best rank-1 tensor formula for minimizing KL divergence (𝑑 = 3 )The best rank-1 approximation of 𝒫 ∈ ℝ𝐼×𝐽×𝐾
is given as>0Normalized vector
depending on only 𝑖

Normalized vector
depending on only 𝑗

By the way,
Frobenius error
minimization
is NP-hard

Normalized vector
depending on only 𝑘which minimizes KL divergence from 𝒫.
We reproduce the result in K.Huang, et al. "Kullback-Leibler principal component for tensors is not NP-hard." ACSSC 2017

A tensor with 𝑑 indices is a joint distribution with 𝑑 random variables.
A vector with only 1 index is an independent distribution with only one random variable.

Rank-1 approximation approximates a joint distribution by a product of independent distributions.
Mean-field approximation : a methodology in physics for reducing a many-body problem to a one-body problem.
38

Mean-field approximation and rank-1 approximation
MFA of Boltzmann-machine1
𝑝 𝒙 =exp ෍ 𝜃𝑖 𝑥𝑖 + ෍ 𝜃𝑖𝑗 𝑥𝑖 𝑥𝑗
𝑍(𝜽)
𝑖

𝑖<𝑗

Bias

Interaction11

𝜂𝑖 = ෍ ⋯ ෍ 𝑥𝑖 𝑝 𝒙
𝑥1 =0

𝑥𝑛=0

𝐷𝐾𝐿 𝑝, 𝑝Ƹ

39

Mean-field approximation and rank-1 approximation
MFA of Boltzmann-machine1
𝑝 𝒙 =exp ෍ 𝜃𝑖 𝑥𝑖 + ෍ 𝜃𝑖𝑗 𝑥𝑖 𝑥𝑗
𝑍(𝜽)

=

𝑖

𝑖<𝑗

Bias

Interaction

1
exp ෍ 𝜃𝑖 𝑥𝑖
𝑍(𝜽)

1

1

𝜂𝑖 = ෍ ⋯ ෍ 𝑥𝑖 𝑝 𝒙
𝑥1 =0

𝑥𝑛=0

= 𝑝 𝑥1 … 𝑝(𝑥𝑛 )𝑖

𝐷𝐾𝐿 𝑝, 𝑝Ƹ

40

Mean-field approximation and rank-1 approximation
MFA of Boltzmann-machine
𝑂 2𝑛
1
𝑝 𝒙 =
exp ෍ 𝜃𝑖 𝑥𝑖 + ෍ 𝜃𝑖𝑗 𝑥𝑖 𝑥𝑗
𝑍(𝜽)
𝑖

𝑖<𝑗

Bias

Interaction

MF equations11

𝜂𝑖 = ෍ ⋯ ෍ 𝑥𝑖 𝑝 𝒙
𝑥1 =0

𝑥𝑛=0

𝐷𝐾𝐿 𝑝Ƹ𝑒 , 𝑝𝜂ҧ𝑖 = sigmoid 𝜃𝑖 + ෍ 𝜃𝑘𝑗 𝜂ҧ𝑘𝑘

𝐷𝐾𝐿 𝑝, 𝑝Ƹ𝑂 2𝑛41

Mean-field approximation and rank-1 approximation
MFA of Boltzmann-machine
𝑂 2
1

1
𝑝 𝒙 =
exp ෍ 𝜃𝑖 𝑥𝑖 + ෍ 𝜃𝑖𝑗 𝑥𝑖 𝑥𝑗𝑍(𝜽)𝑖MF equations1

𝜂𝑖 = ෍ ⋯ ෍ 𝑥𝑖 𝑝 𝒙

𝑖<𝑗

Bias

Rank-1 approximation

𝑛

𝑥1 =0

𝑥𝑛=0

𝑖

𝑗

𝑝𝜃 (𝑖, 𝑗, 𝑘) = exp ෍ ෍ ෍ 𝜃𝑖 ′ 𝑗 ′ 𝑘 ′𝑖 ′ =1 𝑗 ′ =1 𝑘 ′ =1

Interaction

𝐽

𝑘

𝐾

𝜂𝑖11 = ෍ ෍ 𝒫𝑖𝑗 ′ 𝑘 ′𝑗 ′ =1 𝑘 ′ =1

𝐷𝐾𝐿 𝑝Ƹ𝑒 , 𝑝𝜂ҧ𝑖 = sigmoid 𝜃𝑖 + ෍ 𝜃𝑘𝑗 𝜂ҧ𝑘𝑘

𝐷𝐾𝐿 𝑝, 𝑝Ƹ
𝑂 2𝑛

𝐷𝐾𝐿 𝑝, 𝑝Ƹ

Set of products of independent distributions
42

Mean-field approximation and rank-1 approximation
MFA of Boltzmann-machine
𝑂 2
𝑝 𝒙 =

1

1
exp ෍ 𝜃𝑖 𝑥𝑖 + ෍ 𝜃𝑖𝑗 𝑥𝑖 𝑥𝑗𝑍(𝜽)𝑖MF equations1

𝜂𝑖 = ෍ ⋯ ෍ 𝑥𝑖 𝑝 𝒙

𝑖<𝑗

Bias

Rank-1 approximation

𝑛

𝑥1 =0

𝑥𝑛=0

𝐶𝑜𝑚𝑝𝑢𝑡𝑎𝑏𝑙𝑒
𝑖

𝑗

𝐽

𝑘

𝑝𝜃 (𝑖, 𝑗, 𝑘) = exp ෍ ෍ ෍ 𝜃𝑖 ′ 𝑗 ′ 𝑘 ′𝜂𝑖11 = ෍ ෍ 𝒫𝑖𝑗 ′ 𝑘 ′𝑖 ′ =1 𝑗 ′ =1 𝑘 ′ =1

Interaction𝐾𝑗 ′ =1 𝑘 ′ =1

𝐷𝐾𝐿 𝑝Ƹ𝑒 , 𝑝𝜂ҧ𝑖 = sigmoid 𝜃𝑖 + ෍ 𝜃𝑘𝑗 𝜂ҧ𝑘𝑘

𝐷𝐾𝐿 𝑝, 𝑝Ƹ
𝑂 2𝑛

𝐷𝐾𝐿 𝑝, 𝑝Ƹ
𝑂 𝐼𝐽𝐾

Set of products of independent distributions
43

Mean-field approximation and rank-1 approximation

Mean-field
Approximationof BM

Projection onto
e-flat space

Rank-1approximation

Projection onto
e-flat space

Minimizing KL divergence

Minimizing inverse-KL divergence

m-projection

e-projection

impossible𝑛𝑂 2

unique𝜂𝑖 = σ 𝜃𝑖 + ෍ 𝜃𝑘𝑗 𝜂𝑘𝑘not unique

Closed-formula
unique44

Contents□ Motivation, Strategy, and Contributions□ Introduction of log-linear model on DAG□ The best rank-1 approximation formula□ The best rank-1 NMMF

□ Legendre Tucker-Rank Reduction(LTR)□A1GM: faster rank-1 missing NMF

github.com/gkazunii/ Legendre-tucker-rank-reduction

github.com/gkazunii/A1GM□ Theoretical Remarks□ Conclusion
45

12:00

Formulate Tucker rank reduction by relaxing the rank-1 condition
Rank-1 condition (𝜽-representation)
rank 𝒫 = 1 ⟺ its all many−body 𝜃 parameters are 0
𝜃311
𝜃211
𝜃112 𝜃113
𝜃121 𝜃

𝑖𝑗𝑘 = 0

Expand the tensor by focusing on the 𝑚-th axis into a rectangular matrix 𝜃 (𝑚)(mode-𝑚 expansion)

𝜃131

46

Formulate Tucker rank reduction by relaxing the rank-1 condition
Rank-1 condition (𝜽-representation)
rank 𝒫 = 1 ⟺ its all many−body 𝜃 parameters are 0
𝜃311
𝜃211

Expand the tensor by focusing on the 𝑚-th axis into a rectangular matrix 𝜃 (𝑚)(mode-𝑚 expansion)

𝜃112 𝜃113
𝜃121 𝜃

𝑖𝑗𝑘 = 0𝜃131

𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
0

𝜃131
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃113
0
0

0 0
0 0
0 0

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃311
0
0

0
0
0

0
0
0

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
0
0

𝜃121
0
0

0
0
0

0
0
0

𝜃131
0
0

0
0
0

0
0
0

Rank 1,1,1

47

Formulate Tucker rank reduction by relaxing the rank-1 condition
Rank-1 condition (𝜽-representation)
rank 𝒫 = 1 ⟺ its all many−body 𝜃 parameters are 0
𝜃311
𝜃211

Expand the tensor by focusing on the 𝑚-th axis into a rectangular matrix 𝜃 (𝑚)(mode-𝑚 expansion)

𝜃112 𝜃113
𝜃121 𝜃

𝑖𝑗𝑘 = 0𝜃131

𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
0

𝜃131
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃113
0
0

0 0
0 0
0 0

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃311
0
0

0
0
0

0
0
0

Two bingos

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
0
0

𝜃121
0
0

0
0
0

0
0
0

𝜃131
0
0

0
0
0

0
0
0

Two bingos

Rank 1,1,1

Two bingos

48Formulate Tucker rank reduction by relaxing the rank-1 condition
Rank-1 condition (𝜽-representation)
rank 𝒫 = 1 ⟺ its all many−body 𝜃 parameters are 0
𝜃311
𝜃211

Expand the tensor by focusing on the 𝑚-th axis into a rectangular matrix 𝜃 (𝑚)(mode-𝑚 expansion)

𝜃112 𝜃113
𝜃121 𝜃

𝑖𝑗𝑘 = 0The first row and first column are the scaling factors

𝜃131

𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
0

𝜃131
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃113
0
0

0 0
0 0
0 0

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
0
0

𝜃112
0
0

0
0
0

0
0
0

𝜃311
0
0

0
0
0

0
0
0

Two bingos

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
0
0

𝜃121
0
0

0
0
0

0
0
0

𝜃131
0
0

0
0
0

0
0
0

Two bingos

Rank 1,1,1

Two bingos

49

The relationship between bingo and rank
𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
𝜃321

𝜃131
0
𝜃331

𝜃112
0
𝜃312

0
0
𝜃322

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
𝜃321
𝜃331

𝜃112
0
0

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
𝜃312
𝜃313

𝜃121
0
0

0
0

0
0

0
0

𝜃332

𝜃113
0
𝜃313

𝜃323

𝜃333

0
0
0

𝜃312
𝜃322
𝜃332

𝜃311
0
0

0
0
0

𝜃313
𝜃323
𝜃333

No bingo

0
0
0

𝜃321
𝜃322
𝜃323

𝜃131
0
0

0
0
0

𝜃331
𝜃332
𝜃333

No bingo

One bingo

Rank 2,3,3

50

The relationship between bingo and rank
𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
𝜃321

𝜃131
0
𝜃331

𝜃112
0
𝜃312

0
0
𝜃322

0
0
𝜃332

𝜃113
0
𝜃313

0
0

0
0

𝜃323

𝜃333

One bingo

𝜃123
●

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
𝜃321
𝜃331

𝜃112
0
0

0
0
0

𝜃312
𝜃322
𝜃332

𝜃311
0
0

0
0
0

𝜃313
𝜃323
𝜃333

No bingo

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
𝜃312
𝜃313

𝜃121
0
0

0
0
0

𝜃321
𝜃322
𝜃323

𝜃131
0
0

0
0
0

𝜃331
𝜃332
𝜃333

No bingo

Rank 2,3,3

𝒫 Input tensor
𝐷𝐾𝐿 𝒫, 𝒫തm-projection𝒫ത●

Subspace with one bingo in the mode-1 direction ℬ 1

51

The relationship between bingo and rank
𝜃111
𝜃 (1) = 𝜃211
𝜃311

𝜃121
0
𝜃321

𝜃131
0
𝜃331

𝜃112
0
𝜃312

0
0
𝜃322

0
0
𝜃332

𝜃113
0
𝜃313

0
0

0
0

𝜃323

𝜃333

One bingo

𝜃123
●

𝜃111
𝜃 (2) = 𝜃121
𝜃131

𝜃211
0
0

𝜃311
𝜃321
𝜃331

𝜃112
0
0

0
0
0

𝜃312
𝜃322
𝜃332

𝜃311
0
0

0
0
0

𝜃313
𝜃323
𝜃333

No bingo

𝜃111
𝜃 (3) = 𝜃112
𝜃113

𝜃211
0
0

𝜃311
𝜃312
𝜃313

𝜃121
0
0

0
0
0

𝜃321
𝜃322
𝜃323

𝜃131
0
0

0
0
0

𝜃331
𝜃332
𝜃333

No bingo

Rank 2,3,3

𝒫 Input tensor
𝐷𝐾𝐿 𝒫, 𝒫തm-projection𝒫ത●

Subspace with one bingo in the mode-1 direction ℬ 1

Bingo rule (𝑑 = 3 )𝐼 ×𝐼 ×𝐼3

1 2
The mode-𝑘 expansion 𝜃 (𝑘) of the natural parameter of a tensor 𝒫 ∈ ℝ>0⇒ rank 𝒫 ≤ 𝐼1 − 𝑏1 , 𝐼2 − 𝑏2 , 𝐼3 − 𝑏3

has 𝑏𝑘 bingos

52

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.

53

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.

Bingo

Bingo

Bingo

54

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.
STEP2 : Replace the bingo part withthe best rank-1 tensor.

The shaded areas do not change their values in the projection.
55

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.
STEP2 : Replace the bingo part withthe best rank-1 tensor.

Replace the partial tensor in the red box using the best rank-1 approximation formula56

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.
STEP2 : Replace the bingo part withthe best rank-1 tensor.

Replace the partial tensor in the red box using the best rank-1 approximation formula
57

Example: Reduce the rank of (8,8,3) tensor to (5,8,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.
STEP2 : Replace the bingo part withthe best rank-1 tensor.

Replace the partial tensor in the red box using the best rank-1 approximation formula
The best tensor is obtained in the specified bingo space. 😄
There is no guarantee that it is the best rank (5,8,3) approximation.

😢

58

Example: Reduce the rank of (8,8,3) tensor to (5,7,3) or less
𝜃 can be any𝜃 is zero

STEP1 : Choose a bingo location.
STEP2 : Replace the bingo part withthe best rank-1 tensor.

The shaded areas do not change their values in the projection.
59

Experimental results (synthetic data)

LTR is faster with the competitive approximation performance.

60Experimental results (real data)(92, 112, 400)

(9, 9, 512, 512, 3)

LTR is faster with the competitive approximation performance.

61

Contents□ Motivation, Strategy, and Contributions□ Introduction of log-linear model on DAG□ The best rank-1 approximation formula□ The best rank-1 NMMF

□ Legendre Tucker-Rank Reduction(LTR)□A1GM: faster rank-1 missing NMF

github.com/gkazunii/ Legendre-tucker-rank-reduction

github.com/gkazunii/A1GM□ Theoretical Remarks□ Conclusion
62

16:40

Strategy for rank-1 NMF with missing values□ Collect missing values in a corner of matrix to solve as coupled NMF

Missing value

Element-wise product 𝚽 = ቊ0 If 𝐗 𝑖𝑗 is missing𝑖𝑗
1 otherwise

63

Strategy for rank-1 NMF with missing values□ Collect missing values in a corner of matrix to solve as coupled NMF
NMMF (Takeuchi et al., 2013)
EquivalentMissing value

Element-wise product 𝚽 = ቊ0 If 𝐗 𝑖𝑗 is missing𝑖𝑗
1 otherwise64

NMMF, Nonnegative multiple matrix factorization (Takeuchi et al., 2013)
artist

user

user

user

user

user

artisttag

tag

artist65

The best rank-1 approximation of NMMF
artist

user

user

user

user

user

artisttag

tag

artistThe best rank-1 approximation of NMMF

66

Modeling of NMMF

One
To
One

67

One-body and many-body parameters
𝑿, 𝒀, 𝒁 is simultaneously rank-1 decomposable. ⇔ It can be written as 𝒘 ⊗ 𝒉, 𝒂 ⊗ 𝒉, 𝒘 ⊗ 𝒃 .
One-body parameterTwo-body parameter68

Information geometry of rank-1 NMMF
𝑿, 𝒀, 𝒁 is simultaneously rank-1 decomposable. ⇔ It can be written as 𝒘 ⊗ 𝒉, 𝒂 ⊗ 𝒉, 𝒘 ⊗ 𝒃 .
One-body parameterTwo-body parameter

Simultaneous Rank-1 𝜽-condition
Its all two-body 𝜃-parameters are 0.

69

Information geometry of rank-1 NMMF
𝑿, 𝒀, 𝒁 is simultaneously rank-1 decomposable. ⇔ It can be written as 𝒘 ⊗ 𝒉, 𝒂 ⊗ 𝒉, 𝒘 ⊗ 𝒃 .
One-body parameterTwo-body parameter

Simultaneous Rank-1 𝜽-condition
Its all two-body 𝜃-parameters are 0.

Simultaneous Rank-1 𝜼-condition

𝜂𝑖𝑗 = 𝜂𝑖1 𝜂1𝑗
is e-flat. The projection is unique.

70

Find the global optimal solution of rank-1 NMMF
𝑿, 𝒀, 𝒁 is simultaneously rank-1 decomposable. ⇔ It can be written as 𝒘 ⊗ 𝒉, 𝒂 ⊗ 𝒉, 𝒘 ⊗ 𝒃 .
One-body parameterTwo-body parameter

Simultaneous Rank-1 𝜽-condition
Its all two-body 𝜃-parameters are 0.

Simultaneous Rank-1 𝜼-condition

𝜂𝑖𝑗 = 𝜂𝑖1 𝜂1𝑗

The m-projection does not change one-body η-parameterShun-ichi Amari, Information Geometry and Its Applications, 2008, Theorem 11.6

71

Find the global optimal solution of rank-1 NMMF
𝑿, 𝒀, 𝒁 is simultaneously rank-1 decomposable. ⇔ It can be written as 𝒘 ⊗ 𝒉, 𝒂 ⊗ 𝒉, 𝒘 ⊗ 𝒃 .
One-body parameterTwo-body parameter

Simultaneous Rank-1 𝜽-condition
Its all two-body 𝜃-parameters are 0.

Simultaneous Rank-1 𝜼-condition

𝜂𝑖𝑗 = 𝜂𝑖1 𝜂1𝑗

The m-projection does not change one-body η-parameterShun-ichi Amari, Information Geometry and Its Applications, 2008, Theorem 11.6

All 𝜼-parameters after the projection are identified.

72

19:20

Rank-1 NMF with missing values□ NMMF can be viewed as a special case of NMF with missing values.
Equivalent73

Rank-1 NMF with missing values□ NMMF can be viewed as a special case of NMF with missing values.
Equivalent□ NMF is homogeneous for row and column permutations

74

A1GM: Algorithm

Step 1 : Gather missing values in the bottom right.

Step 2 : Use the formula of the best rank-1 NMMF.
Step 3 : Repermutate

Find exact solution 🤔❓

75

Examples that permutations cannot collect missing values into corners

76Add missing values to solve the problem as NMMF

77

Add missing values to solve the problem as NMMF

Reconstruction error worsens 😢
78

Add missing values to solve the problem as NMMF

Reconstruction error worsens 😢

Gain in efficiency 😀

79🙆Data that A1GM is good at and not good at🙅🙅 Missing values are evenly distributed in each row and column.

80🙆Data that A1GM is good at and not good at🙅🙅 Missing values are evenly distributed in each row and column.

🙆 Missing are heavily distributed in certain rows and columns.
Missing values tend to be in certain columns in some real datasets.
ex) disconnected sensing device,
optional answer field in questionnaire form

81

A1GM: Algorithm

Step 1 : Increase the number of missing values.
Step 2 : Gather missing values in the bottom right.

Step 3 : Use the formula of rank-1 NMMF and repermutate.

82

Experiments on real data
□ A1GM is compared with gradient-based KL-WNMF- Relative runtime < 1 means A1GM is faster than KL-WNMF.
- Relative error > 1 means worse reconstruction error of A1GM than KL-WNMF.
- Increase rate is the ratio of # missing values after addition of missing values at step1.

×5 – 10times faster!

Find
the best solution

Add missing values.
Accuracy decreases.

83

Contents□ Motivation, Strategy, and Contributions□ Introduction of log-linear model on DAG□ The best rank-1 approximation formula□ The best rank-1 NMMF

□ Legendre Tucker-Rank Reduction(LTR)□A1GM: faster rank-1 missing NMF

github.com/gkazunii/Legendre-tucker-rank-reduction

github.com/gkazunii/A1GM□ Theoretical Remarks□ Conclusion
84

22:30

Theoretical Remarks 1 : Extended NMMF.
□ The rank of weight matrix is 2 after adding missing values.

𝐗𝚽rank 𝚽 = 2

𝐗

0 If 𝐗 𝑖𝑗 is missing
𝚽𝑖𝑗 = ቊ1 otherwise𝚽rank 𝚽 = 2

85

Theoretical Remarks 1 : Extended NMMF.
□ The rank of weight matrix is 2 after adding missing values.

𝐗𝚽rank 𝚽 = 2

𝐗

0 If 𝐗 𝑖𝑗 is missing
𝚽𝑖𝑗 = ቊ1 otherwise𝚽rank 𝚽 = 2□ Can we exactly solve rank-1 NMF if the rank(Φ) = 2?

rank 𝚽 = 2

86

Theoretical Remarks 1 : Extended NMMF.
The best rank-1 approximation of extended NMMF

87

Theoretical Remarks 1 : Extended NMMF.
The best rank-1 approximation of extended NMMF

Equivalent𝚽𝑖𝑗 = ቊ0 If 𝐗 𝑖𝑗 is missing
1 otherwise
88

Theoretical Remarks 1 : Extended NMMF.
The best rank-1 approximation of extended NMMF

Equivalent
Permutation

𝚽𝑖𝑗 = ቊIf rank(𝚽) ≦2,
the matrix can be transformed into the form

0 If 𝐗 𝑖𝑗 is missing
1 otherwise

We can exactly solve rank-1 NMF with missing values by permutation if rank(𝚽)≦2.
89Theoretical Remarks 2 : Connection to balancing.
□ Matrix Balancing

Mahito Sugiyama, Hiroyuki Nakahara and Koji Tsuda"Tensor balancing on statistical manifold“(2017) ICML.

Balanced matrix(Doubly stochastic matrix)

Transform

90

Theoretical Remarks 2 : Connection to balancing.
□ Matrix Balancing

Mahito Sugiyama, Hiroyuki Nakahara and Koji Tsuda"Tensor balancing on statistical manifold“(2017) ICML.

Balancing condition

Balanced matrix(Doubly stochastic matrix)

Transform

91

Theoretical Remarks 2 : Connection to balancing.
□ Matrix Balancing

Balancing condition

Balanced matrix(Doubly stochastic matrix)

Transform

Rank-1 condition
Its all many-body𝜃-parameters are 0.

Balanced rank-1 matrix is unique.
92

Conclusion
Data structure

DAG

Infor-Geo□ Closed Formula of the Best Rank-1 NMMF

The best rank-1 approximation for NMMF

□ Describe low-rank condition using (𝜃,𝜂)Rank-1 condition (𝜼-representation)
𝜂ҧ 𝑖𝑗𝑘 = 𝜂ҧ 𝑖11 𝜂ҧ1𝑗1 𝜂ҧ11𝑘
Rank-1 condition (𝜽-representation)
ҧ are 0All many body 𝜃𝑖𝑗𝑘

□ A1GM: Faster Rank-1 NMF with missing values

9393