JP5684055B2 - Time-related data analysis apparatus, method, and program - Google Patents

Description

  The present invention relates to a time relationship data analysis device, method, and program, and more particularly, to a time relationship data analysis device, method, and program for analyzing a time series of relationship data between objects.

  Currently, a large amount of data is available through the Internet, and information analysis and business using such a large-scale data set are actively performed. However, it is difficult for humans to grasp the contents of such a large amount of data.

  In particular, recently, attention has been focused on “relation data” indicating relationships between objects such as human relationships and business relationships between companies. By using such relationship data, for example, it is possible to realize purchase item recommendation with higher accuracy from the relationship between friends.

  Many methods for analyzing relational data between objects have been proposed so far. Many methods focus on the fact that the relationship between the two can be expressed in the form of a matrix, and structure the matrix. Non-Patent Document 1 and Non-Patent Document 2 perform clustering of objects having similar connection destinations when information on the presence / absence of a relationship between two objects is given. In addition, since there are many people and items having similar tendencies (patterns) in these relation data, it is known that the rank is low when the relation data is viewed as a matrix. Non-Patent Document 6 is a method for estimating an unobserved part of a partially observed matrix using this property as an assumption.

  The problem with these methods is that, for example, relationships between friends and transactions between companies change over time. Therefore, in order to realize more accurate relational data analysis, a technique for observing a relation between certain objects at a plurality of times and analyzing temporal relational data including the change is required. For this purpose, a method for analyzing data in which a matrix representing a relationship is defined for each time is required, not just analysis of matrix data. Many methods for analyzing time-related data have not yet been proposed. The methods described in Non-Patent Literature 3 and Non-Patent Literature 4 are methods that have been extended so that Non-Patent Literature 1 and Non-Patent Literature 2 can be applied when time-varying. In the method described in Non-Patent Document 3, modeling is performed assuming a situation in which the relationship changes continuously by using a state space model. On the other hand, Non-Patent Document 4 proposes a time evolution model that can cope with a rapid time change by using an HMM. The method described in Non-Patent Document 5 models time-related data as a tensor and models time changes using a kind of autoregressive model.

K. Nowicki and T. A. B. Snijders. Estimation and prediction for stochastic blockstructures. Journal of the American Statistical Association, 96 (455): 1077-1087, 2001. C. Kemp, J. B. Tenenbaum, T. L. Griffiths, T. Yamada, and N. Ueda. Learning systems of concepts with an infinite relational model.In Proceedings of the 21st National Conference on Artificial Intelligence (AAAI), 2006. T. Yang, Y. Chi, S. Zhu, Y. Gong, and R. Jin.A Bayesian approach toward finding communities and their evolutions in dynamic social networks.In Proceedings of SIAM International Conference on Data Mining (SDM), 2009. K. Ishiguro, T. Iwata, N. Ueda and J. Tenenbaum. Dynamic Infinite Relational Model for Time-varying Relational Data Analysis. Advances in Neural Information Processing Systems 23, 2010. L. Xiong, X. Chen, T.-K. Huang, J. Schneider, and J. G. Carbonell.Temporal Collaborative Filtering with Bayesian Probabilistic Tensor Factorization.In Proceedings of SIAM International Conference on Data Mining (SDM), 2010. N. Srebro, J. D. M. Rennie, and T. S. Jaakkola. Maximum-Margin Matrix Factorization. In Advances in Neural Information Processing Systems 17, 2005.

  Each of the above methods has the following problems.

  First, the methods described in Non-Patent Document 1, Non-Patent Document 2, and Non-Patent Document 6 have a problem that it is impossible to analyze time-related data.

  The methods described in Non-Patent Document 3 and Non-Patent Document 5 can divide time-related data into patterns, but have a problem that the number of patterns needs to be determined in advance.

  The method described in Non-Patent Document 4 is a time-related data analysis method that can automatically determine the number of patterns. However, only the presence or absence of the relationship between objects is the object of analysis. Therefore, it is unsuitable for analyzing a more accurate relationship, such as a change in the strength of the relationship. Moreover, it is impossible to estimate the time when no observation is made by this method. For example, when time-related data such as an input-output table surveyed every five years is given, it is impossible in principle to restore relationships in other periods (such as the third year) and predict future relationships. There is a problem that there is.

  The present invention has been made in order to solve the above-described problem, and analyzes the relationship data between objects that change over time to obtain the strength of the relationship between unobserved objects. It is an object to provide a method and a program.

In order to achieve the above object, the time relationship data analysis device represents a time series of relationship data indicating the strength of the relationship between each of the plurality of first objects and each of the plurality of second objects, and A time relationship data analysis apparatus that analyzes time relationship data expressed by a tensor and obtains an unknown tensor X representing the time relationship data, wherein the relationship between the first object and the second object is strong. Represents the strength of regularization relative to each mode k (k = 1, 2, 3) of the observation vector y consisting of M observation values, the regularization constant λ which is a positive real number, and the tensor X The storage means storing the regularization balance constant γ _k , and the observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) are acquired from the storage means, and Represented by the formula (I) Analyzing means for obtaining the tensor X that minimizes the objective function.

Here, n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the time-related data, and x is an element of the tensor X. Ω (X) is an operation for extracting M elements corresponding to the observation value of the observation vector from the tensor X. mat ₁ is a mode 1 matrix that returns n ₁ × n ₂ n ₃ , mat ₂ is a mode 2 matrix that returns n ₂ × n ₃ n ₁ , and mat ₃ is n ₃ × n ₁ n Mode 3 matrixing that returns ₂ . ‖‖ _tr represents the linear sum of the singular values of the matrix.

The temporal relationship data analysis method according to the present invention includes an observation vector y composed of M observation values of the strength of the relationship between the first object and the second object, a regularization constant λ that is a positive real number, and a tensor. In a time-related data analysis apparatus including a storage unit that stores a regularization balance constant γ _k representing a relative regularization strength relative to each mode k (k = 1, 2, 3) of X. A time series of relation data indicating the strength of the relation between each of the objects and each of the plurality of second objects is represented, and the time relation data expressed by a tensor is analyzed to represent the time relation data a temporal relationship data analysis method for determining the unknown tensor X, by the analyzing means, from the storage unit, the observation vector y, the regularization constant lambda, the regularization balance constant γ _{k (k} = 1,2 3) Get the is characterized by determining said tensor X that minimizes the objective function represented by the above formula (I).

According to the present invention, the analysis unit obtains the observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) from the storage unit, I) Find the tensor X that minimizes the objective function expressed by the equation.

  In this way, by analyzing the relationship data between the objects that change over time and obtaining the unknown tensor X representing the time relationship data so as to minimize the objective function, the strength of the relationship between the unobserved objects is obtained. Can be requested.

The storage means according to the present invention includes an auxiliary variable Z _k (k = 1, 2, 3), a value x of each element of the tensor X, a Lagrange multiplier A _k (k = 1, 2, 3), a first object A pattern base U _k (k = 1, 2, 3) in each of a direction, a second object direction, and a time direction, the first object changing according to the pattern base U _k (k = 1, 2, 3), A relationship pattern V _k (k = 1, 2, 3) representing the strength of the relationship between the two objects and each pair of time is further stored, and the analysis means includes an auxiliary coefficient Zk and a Lagrange multiplier for each mode k. Based on A _k , the predicted value of the tensor ^ X is calculated according to the following equation (II), and the tensor is calculated according to the following equation (III) based on the calculated predicted value of the tensor ^ X and the observed vector y: The value x of each element of X is updated and stored in the storage means. A main variable updating means for, based on the Lagrange multipliers A _k of the tensor X and modes k, for each mode k, in accordance with the following formula (IV), to calculate the predicted value ^ Z _k of auxiliary variables, calculates On the basis of the predicted value ^ Z _k of the auxiliary variable, for each mode k, singular value decomposition is performed according to the following equation (V) to obtain a diagonal matrix S _{k in} which singular values are arranged, The pattern base U _k and the relation pattern V _k are updated and stored in the storage means. Based on the obtained diagonal matrix S _k , the pattern base U _k , and the relation pattern V _k , the following According to the formula (VI), the auxiliary variable Z _k is updated and stored in the storage means, and on the basis of the tensor X and auxiliary variable Z _k , for each mode k, the following (VII) According to the formula A Lagrange multiplier updating means for updating the serial Lagrange multiplier A _k stored in the storage means, updated by the main variable updating means updating by the auxiliary variable updating means, by repeating the updating by the Lagrangian multiplier update means, a predetermined Analysis end condition determining means for outputting the tensor X stored in the storage means when the analysis end condition is satisfied.

  However, ￣Ω (X) is an operation for extracting an element other than the observed value from the tensor X, and η is a constant.

  The temporal relationship data analysis apparatus according to the present invention predicts the relationship between the first object i and the second object j according to the following equation (VIII) based on the tensor X obtained by the analysis means. A relational data predicting means for calculating a predicted value X ′ (i, j) of the relational data at the target time t ′ can be further included.

  However, the times when the observed values are observed before and after the prediction target time t ′ are T (k) and T (k + 1).

  The program according to the present invention is a program for causing a computer to function as each means of the time-related data analysis apparatus.

  As described above, according to the time-related data analysis apparatus, method, and program of the present invention, the time-related data is represented so as to minimize the objective function by analyzing the relationship data between the objects that change with time. By obtaining the unknown tensor X, it is possible to obtain the effect that the strength of the relationship between unobserved objects can be obtained.

It is a schematic diagram regarding the analysis of time-related data. It is the schematic which shows the structure of the time-related data analysis apparatus which concerns on embodiment of this invention. It is a figure for demonstrating mode 1 matrixing and mode 1 tensorization. It is a figure for demonstrating mode 2 matrixing and mode 2 tensorization. It is a figure for demonstrating mode 3 matrix formation and mode 3 tensorization. It is a graph which shows the singular value obtained by singular value decomposition. It is a flowchart which shows the content of the time relation data analysis processing routine in the time relation data analysis apparatus which concerns on embodiment of this invention. It is a flowchart which shows the content of the process which updates a main variable. It is a flowchart which shows the content of the process which updates an auxiliary variable, a pattern base, and a relationship pattern. It is a flowchart which shows the content of the process which updates a Lagrange multiplier. It is a flowchart which shows the content of the process which calculates the error of KKT conditions. It is a graph which shows the singular value obtained by experiment. It is a figure which shows the inverse matrix coefficient in 1995 and the estimated inverse matrix coefficient obtained by experiment. It is a figure which shows the singular vector and singular value of the 1st component accompanying the singular value decomposition regarding the mode 3 obtained by experiment. It is a figure which shows the singular vector and singular value of the 2nd component accompanying the singular value decomposition regarding the mode 3 obtained by experiment.

  Hereinafter, embodiments of the present invention will be described in detail with reference to the drawings.

<Overview>
First, time-related data used in the time-related data analysis apparatus proposed in the present invention will be described. FIG. 1 is a diagram schematically showing the concept of the proposed time-related data analysis apparatus.

  In FIG. 1, on the left side, the relationship between the user (purchaser) and the product (whether or not it is bought) is expressed as a tensor along the time axis. The black cube represents the observed relationship, and the other elements of the tensor are unknown. The right side of FIG. 1 shows a low rank decomposition of the left partial observation tensor. Three thin rectangular parallelepipeds surrounding the center represent products, users, and temporal changes with a small number of patterns, respectively, and a rectangular parallelepiped with a small center represents an interaction between patterns. By decomposing in this way, the unknown part of the left tensor can be estimated.

<System configuration>
Next, an embodiment of the present invention will be described by taking as an example the case where the present invention is applied to a time-related data analysis apparatus that obtains and outputs unknown tensors representing time-related data with observation data as an input.

  As shown in FIG. 2, the time-related data analysis device according to the present embodiment includes an input unit 1 that receives input of observation values of time-related data, an analysis unit 2 that analyzes the observed time-related data, A storage unit 3 that stores variables and constants, a predicted value calculation unit 4, and an output unit 5 that outputs an analysis result are provided.

  The input unit 1 is realized by an input device such as a known video camera, microphone, or storage device. The output unit 5 is implemented by a display, a printer, a magnetic disk, or the like.

  The analysis unit 2, the storage unit 3, and the predicted value calculation unit 4 include a CPU (Central Processing Unit), a RAM (Random Access Memory), and a ROM that stores a program for executing a time-related data analysis processing routine described later. (Read Only Memory).

  The analysis unit 2 includes an initial setting unit 21, a variable calculation unit 22, and an analysis end condition determination unit 27. The variable calculation unit 22 includes a main variable update unit 23, an auxiliary variable update unit 24, a Lagrange multiplier update unit 25, and an error calculation unit 26.

    The initial setting unit 21 sets a value for each element of the storage unit 3 from the input observation data.

    The variable calculation unit 22 updates each variable using information stored in the storage unit 3.

    The main variable update unit 23 updates the main variable using information stored in the storage unit 3.

  The auxiliary variable updating unit 24 uses the information stored in the storage unit 3 to update the auxiliary variable, the pattern base, and the relation pattern.

  The Lagrange multiplier update unit 25 uses the information stored in the storage unit 3 to update the Lagrange multiplier.

    The error calculation unit 26 calculates an error using the information stored in the storage unit 3.

    The analysis end condition determination unit 27 determines whether or not the analysis can be ended using the numerical values calculated using the variable calculation unit 22 and stored in the storage unit 3.

  The storage unit 3 includes a constant 31, an observation value 32, an index set 33, a main variable 34, an auxiliary variable 35, a Lagrange multiplier 36, an error 37, a pattern base 38, and a relation pattern 39.

  The predicted value calculation unit 4 performs prediction of related data at the specified predicted time t ′.

<Data structure>
n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the assumed time-related data. The relationship data represents a connection relationship between two types of objects. In the present embodiment, one of the two types of objects constituting the relationship is referred to as a starting object and the other is referred to as an end object. For example, in the example of FIG. 1, the product is a start object and the user is an end object. Further, in the present embodiment, n ₁ and n ₂ are the number of objects at the start and end of the relationship, respectively, and n ₃ is the total number of steps at the time when the observation occurred. Further, N is the number of all possible elements, and N = n ₁ × n ₂ × n ₃ . The axis of the object that is the starting point is called the first dimension of the tensor, the axis of the object that is the end point is called the second dimension of the tensor, and the axis of time is called the third dimension of the tensor.

T is a vector that stores the times of n ₃ observed time points. Assume that the elements of T are arranged in ascending order.

The observation vector y is observed relational data and is expressed by a vertical vector of length M. M is the number of elements of y, _{y m} represents the observed value of the m-th y (1 ≦ m ≦ M) . Note that the observed values y _m in the present embodiment allows for all of the real numbers. Therefore, it can be understood that this value represents the strength of a certain relationship. On the other hand, in the non-patent document 3 and the non-patent document 4 described above, since the observed relational data only allows binary values of 1 or 0, only the presence or absence of the relation can be expressed, and the analysis considering the strength of the relation is impossible. It is.

I, J, and K are each an index set 33 of length M, and (i _m , j _m , k _m ) in which the m-th elements are arranged indicates the m-th observed element of the unknown tensor. That m th observation value is the relationship with i _m-th origin object and the j _m-th time k _m endpoint object.

Similarly, I ^c , J ^c , and K ^c are index sets 33 each having a length (N−M) in which m-th elements are arranged (i ^c _m , j ^c _m , k ^c _m ). Indicates the mth unobserved element of the unknown tensor.

  λ is a positive real number called a regularization constant. Represents the magnitude of the regularization term effect in decomposition / analysis calculations.

  The M-dimensional observation vector y is stored in the observation value 32.

The main variable 34 stores an unknown tensor X composed of n ₁ × n ₂ × n ₃ elements. An appropriate value is input to each element of X as an initial value. In general, all values should be 0. The unknown tensor X corresponds to the tensor shown in FIG. 1, and finally, the element (i, j, k) is between the start object i and the end object j observed at time k. A value indicating the strength of the relationship is stored (overwritten).

The following four operations are defined: observation element extraction Ω, unobserved element extraction ΩΩ, observation element restoration Ω ^T , mode k matrixization mat _k , and mode k tensorization tensor _k . These operations are expressed as matrix operations in the following mathematical formula, but actually do not need to be held as a matrix. First, the observation element extraction Ω takes a tensor consisting of n ₁ × n ₂ × n ₃ elements as an argument and returns an M-dimensional vector. The m-th element of the return value is the first of the tensors that are arguments. (i _m , j _m , k _m ) elements. Similarly, the unobserved element extraction ￣Ω takes a tensor consisting of n ₁ × n ₂ × n ₃ elements as an argument and returns a (N−M) -dimensional vector. The m-th element of the return value is , Is the (i ^c _m , j ^c _m , k ^c _m ) element of the argument tensor. The observation element restoration Ω ^T is defined as the inverse operation of the observation element extraction Ω, and returns a tensor consisting of n ₁ × n ₂ × n ₃ elements with an M-dimensional vector as an argument. The m th element of the argument Is stored in the (i _m , j _m , k _m ) element of the tensor that is the return value. Return elements that do not correspond to observations are all zero. There are three types of mode k matrix mat- _k , k = 1, 2, and 3, each of which has an n ₁ × n ₂ × n ₃ element tensor as an argument, an n ₁ × n ₂ n ₃ matrix, n Returns a ₂ × n ₃ n ₁ matrix and an n ₃ × n ₁ n ₂ matrix. Specifically, the mode 1 matrixing mat ₁ generates the first mode of the given tensor as shown in FIG. 3A, that is, the vector is extracted so that the starting objects are vertically continuous and arranged in a matrix. It is an operation to do. The first mode tensorization tensor- ₁ is the reverse operation. Similarly, the second mode matrixing mat- ₂ is an operation for taking out a vector so that the end point objects are vertically continuous as shown in FIG. 3B and generating a matrix arranged horizontally, and the second mode tensor. The tensor ₂ is the reverse operation. As shown in FIG. 3C, the third mode matrixing mat ₃ is an operation for extracting vectors such that time points are vertically continuous and generating a matrix arranged side by side, and the third mode tensorization tensor- ₃ is The reverse operation.

<Input unit 1>
The input unit 1 accepts input of numerical data required in advance for analysis. The input data is as follows.

n ₁ , n ₂ , n ₃ , N, M, I, J, K, I ^c , J ^c , K ^c , y, λ, γ ₁ , γ ₂ , γ ₃ , η, L, E, t ′

n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the assumed time-related data. Let N be all possible elements and N = n ₁ × n ₂ × n ₃ .

y is the observed relational data, and is an observation vector expressed by a vertical vector of length M. M is the number of elements of y, y _m denotes the m th observation.

I, J, and K are indices of observation values, and represent the relationship between M observation values y and unknown tensors X. Similarly, I ^c , J ^c , and K ^c that are indices of (N−M) unobserved values are also input.

λ is a positive real number called a regularization constant. Represents the magnitude of the regularization term effect in decomposition / analysis calculations. γ ₁ , γ ₂ , and γ ₃ are positive constants representing the relative regularization strength of each mode k of the unknown tensor X. Here, these are called regularization balance constants. In addition, the solution obtained by multiplying the regularization constant λ by a certain constant and dividing the regularization balance constants γ ₁ , γ ₂ , γ ₃ by the constants is the same. Therefore, for normalization, γ ₁ + γ ₂ + γ ₃ = 1.

η is a step size, which is a positive real number. It is a constant that determines the amount of update of the Lagrange multipliers A ₁ , A ₂ , A _{3 in} the decomposition / analysis calculation, and represents the strength of the penalty for the difference between the main variable and the auxiliary variable.

  L is a natural number representing the maximum number of repetitions in the variable calculation unit 22. The number of repetitions is represented by l.

  E is a constant representing the upper limit of the allowable error in the analysis end condition determination unit 27.

  t ′ is a real number representing the index of the time to be predicted. Note that t 'may not be input when prediction is unnecessary.

<Analysis unit 2>
The analysis unit 2 obtains an unknown tensor X that is a main variable based on the input observation values and constants. The analysis unit 2 determines the main variable X by approximately minimizing an objective function expressed by the following equation (1).

Here, the first term of the equation (1) is a loss term, and the square of the difference between the observation vector y and the M-dimensional vector Ω (X) obtained by extracting the element corresponding to the observation value from the unknown tensor X Defined as sum. However, when there is a missing value in the observed value, the unknown tensor X cannot be uniquely determined only by the first term. The second term is called the regularization term and expresses the assumption that the rank of the unknown tensor is low. Here, || Z _k || _tr is called a trace norm and is defined as a linear sum of singular values of the matrix Z _k . The fact that the unknown tensor is decomposed into a small number of patterns and the interactions between them as shown in FIG. 1 corresponds to the fact that the matrix obtained by the matrixing operation mat- _k has a low rank. Note (for example, see Non-Patent Document 8 (TG Kolda and BW Bader. Tensor decompositions and applications. SIAM Review, 51 (3): 455-500, 2009.)). To explain why trace norm regularization corresponds to the hypothesis of low rank, only regularization terms corresponding to the first mode are considered for simplicity when all elements are observed (ie, γ ₁ = 1, γ ₂ = 0, γ ₃ = 0). At this time, the above optimization problem can be rewritten as the following equation (2).

Here, X and Y are assumed to be matrices obtained by mode 1 matrixization of the unknown tensor X and the observation vector y, respectively. It is known that the solution of the above optimization problem can be analytically written as the following equation (3) where the singular value decomposition of the matrix Y is Y = USV ^T (for example, Non-Patent Document 9 (R Tomioka, T. Suzuki, and M. Sugiyama. Augmented lagrangian methods for learning, selecting, and combining features. In S. Sra, S. Nowozin, and SJ Wright, editors, Optimization for Machine Learning. MIT Press, 2011.) See).

  Here, S is a diagonal matrix in which singular values of the matrix Y are arranged. U and V are matrices in which each column is a singular vector corresponding to a singular value. Further, the max operation is applied to each element of the matrix. When this optimization problem is solved from the above equation (3), as shown in FIG. 4, the singular values below the regularization constant λ among the singular values of the given matrix Y are discarded as noise, and the signal is sufficiently larger than the noise. It turns out that an ingredient can be taken out. Further, the above operation is invariant even when an orthogonal matrix is applied from both sides of the observed value Y, and is robust to these operations.

  For example, Non-Patent Document 7 (S. Boyd, N. Parikh, E. Chu, B. Peleato, and J. Eckstein. Distributed optimization and statistical learning via the alternating direction method of multipliers. It can be solved using an Alternative Direction Method of Multipliers (ADMM) method such as www.stanford.edu/~boyd/papers/distr_opt_stat_learning_admm.html). Furthermore, since the above optimization problem is a convex optimization problem (refer to Non-Patent Document 10 (S. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press. 2004.) for general theory on convex optimization), the ADMM method is If it converges, it is the global optimal solution of the above optimization problem. The convergence of the ADMM method is determined using the error of the KKT condition as described in detail below. Thus, it is not necessary to specify how many patterns (dimensions) the tensor is decomposed into, and these can be automatically determined through optimization.

<Initial setting unit 21>
The initial setting unit 21 performs initial setting of each element stored in the storage unit 3 based on information input from the input unit 1.

First, n ₁ , n ₂ , n ₃ , N, M, y, λ, γ ₁ , γ ₂ , γ ₃ , η, L, E, and t ′ are stored in the constant 31 from the input unit 1. Also, the number of repetitions l is initialized to 1.

The M-dimensional observation vector y is stored in the observation value 32, and the index I, J, K of the observation value and the index I ^c , J ^c , K ^c of the unobserved value are stored in the index set 33.

The main variable 34 stores an unknown tensor X composed of n ₁ × n ₂ × n ₃ elements. An appropriate value is input to each element of X as an initial value. In general, all values should be 0.

-Z ₁ , Z ₂ , Z ₃ are stored in the auxiliary variable 35 and expressed as a matrix of size n ₁ × n ₂ n ₃ , n ₂ × n ₃ n _1, n ₃ × n ₁ n ₂ , respectively. . The vectorization is written as z ₁ , z ₂ , and z ₃ , respectively. The auxiliary variable is obtained by relaxing restrictions on the main variable x in the variable calculation unit 22. Specifically, Z _k is an amount that coincides with the mode k matrixing of the tensor X as optimization proceeds. Appropriate values are input to the auxiliary variables Z ₁ , Z ₂ , and Z ₃ as initial values. In general, all values should be 0.

The Lagrange multiplier 36, Lagrange multiplier A _1, A _2, A ₃ are stored, Lagrange multiplier A _1, A _2, A ₃ is a matrix with the same dimensions as the auxiliary variables _{Z 1, Z 2, Z 3} . The Lagrange multiplier has the role of making the main variable X coincide with the auxiliary variables Z ₁ , Z ₂ , and Z ₃ . As initial values, appropriate values are input to Lagrange multipliers A ₁ , A ₂ , and A ₃ . In general, all values should be 0.

The error _{37, ε 1, ε 2, ···} , ε T is _{stored, ε 1, ε 2, ···} , the epsilon _T represents the error calculated at the error calculating unit 26.

In the pattern base 38, pattern bases U ₁ , U ₂ , U ₃ are stored. Pattern bases U ₁ , U ₂ , and U ₃ are matrices of n ₁ × r ₁ , n ₂ × r ₂ , and n ₃ × r ₃ , respectively, and the origin object direction, end object direction, And express the pattern in the time direction. In the relationship pattern 39, relationship patterns V ₁ , V ₂ , and V ₃ are stored. The relationship patterns V ₁ , V ₂ , V ₃ are matrices of n ₂ n ₃ × r ₁ , n ₃ n ₁ × r ₂ , n ₁ n ₂ × r ₃ , respectively, and are endpoint objects that change according to the pattern base Uk and It represents the strength of the relationship between time, the strength of the relationship between the origin object and time, and the strength of the relationship between the origin object and the end object. r ₁ , r ₂ , and r ₃ each represent the number of elements of the pattern. This value is automatically determined as the tensor decomposition analysis proceeds. As initial values, appropriate values are input to the pattern bases U ₁ , U ₂ , U ₃ , the relationship patterns V ₁ , V ₂ , V ₃ , and the number of pattern elements r ₁ , r ₂ , r ₃ . In general, all values of the pattern base should be zero. For r ₁ , r ₂ , and r _3, an appropriate natural number is selected.

<Variable calculation unit 22>
The variable calculation unit 22 estimates the main variable based on the ADMM method using each element stored in the storage unit 3. In the actual calculation, the analysis processing for performing each processing of the main variable update unit 23, the auxiliary variable update unit 24, the Lagrange multiplier update unit 25, and the error calculation unit 26 in this order is repeatedly performed.

  The analysis end condition determination unit 27 determines the end of the analysis every repetition, and if the analysis is not ended, the number of repetitions t is incremented by 1 and the analysis process is performed again.

<Main variable update unit 23>
The main variable updating unit 23 updates the main variable X using the following equations (4) and (5).

First, the main variable updating unit 23 constructs a predicted value ^ X of the main variable using the above equation (4). The predicted value ^ X is a tensor having the same size as the main variable, and is the average of the tensors _obtained by subtracting the Lagrange multipliers A _k from the auxiliary variables Z _k for each mode k.

  Next, the main variable updating unit 23 updates the main variable X as shown in the above equation (5) using the configured predicted value ^ X of the main variable. In the above equation (5), for the element for which the observed value is obtained, the new main variable value is updated so as to be equal to the predicted value and the observed value divided internally by a ratio of 3λη: 1. Yes. On the other hand, for elements for which no observed value is obtained, the value of the new main variable is updated so as to be the predicted value itself. The obtained main variable X is stored in the main variable 34 of the storage unit 3, and the processing of the main variable update unit 23 is completed.

<Auxiliary variable update unit 24>
The auxiliary variable update unit 24 uses the following equations (6), (7), and (8) to calculate auxiliary variables Z ₁ , Z ₂ , Z ₃ , pattern bases U ₁ , U ₂ , U ₃ , And the relationship patterns V ₁ , V ₂ , V ₃ are updated.

First, as shown in the above equation (6), the predicted value ^ Z _k of the auxiliary variable is configured using the amount stored in the storage unit 3. The auxiliary variable predictive value _{{circumflex over (} Z) _{} k} is obtained by adding a Lagrange multiplier Ak to a matrix of the main variable X for each mode k.

Next, the singular value decomposition is performed on the predicted value _{{circumflex over} (Z) _{} k} of the auxiliary variable configured as described above according to the equation (7). When singular value decomposition is performed, U _k , S _k , and V _k matrices are obtained for the k-th predicted value ^ Z _k . Of these, S _k is a diagonal matrix in which singular values are arranged. U _k and V _k are singular vectors in which each column corresponds to a singular value. However, as the singular value and the singular vector, only the singular value having a value of γ _k / η or more and the corresponding singular vector need be calculated for the operation of the following equation (8).

Finally, the auxiliary variable Z _k is updated according to the above equation (8) using the singular value and the singular vector obtained by the above equation (7). Further, the new pattern base U _k is obtained by excluding the column that becomes 0 by the max function. That is, the number of elements r _k of the pattern may vary in each iteration by the results of a singular value decomposition.

The obtained auxiliary variable Z _k is vectorized and stored in the auxiliary variable 35 of the storage unit 3, the pattern base U _k is stored in the pattern base 38 of the storage unit 3, and the relation pattern V _k is further stored in the relationship of the storage unit 3. The data is stored in the pattern 39, and the analysis process of the auxiliary variable update unit 24 is finished.

<Lagrange multiplier update unit 25>
The Lagrange multiplier updating unit 25 updates the Lagrange multipliers A ₁ , A ₂ , and A ₃ using the following equation (9).

In the above equation (9), the Lagrange multipliers A ₁ , A ₂ , A ₃ are updated by the amount of divergence between the matrix of the main variable X and the auxiliary variables Z ₁ , Z ₂ , Z ₃ . The obtained Lagrangian multiplier is stored in the Lagrange multiplier 36 of the storage unit 3, and the analysis process of the Lagrange multiplier update unit 25 is completed.

<Error calculation unit 26>
The error calculator 26 calculates a KKT condition error. The KKT conditions are the same as those described in Non-Patent Document 10, for example, and will not be described.

The error ε _{t of the} KKT condition is calculated as the maximum value of the _two error terms ε ₁ and ε ₂ as shown in the following equation (10).

However, when the regularization constant λ is a positive value, the error term ε ₁ is calculated according to the following equation (11). On the other hand, when the regularization constant λ is zero, the error term ε ₁ is calculated according to the following equation (12).

The error term ε ₂ is calculated according to the following equation (13).

  The obtained error ε is stored in a location corresponding to the error 37 of the storage unit 3, and the analysis process of the error calculation unit 26 is completed.

<Analysis End Condition Determination Unit 27>
The update of the main variable, auxiliary variable, and Lagrange multiplier is performed iteratively. After each update, the main variable, auxiliary variable, and Lagrange multiplier are updated, and then the analysis end condition determination unit 27 performs the KKT condition. Is used to determine if the iteration should be terminated. Specifically, if the error in the KKT condition is equal to or smaller than a predetermined constant ε or the number of iterations reaches a predetermined number L, the iteration is terminated assuming that the analysis termination condition is satisfied. When continuing the iteration, the number of iterations t is incremented by 1, and the processing returns to the main variable updating unit 23.

  In the element (i, j, k) of the unknown tensor X when the iteration is completed by the analysis end condition determination unit 27, the strength of the relationship between the start object i and the end object j observed at time k is set. The value to represent will be stored. Thereby, it is possible to acquire data (unknown tensor X) obtained by restoring relational data between objects that are affected by noise or that cannot be observed. In addition, by paying attention to each axis (dimension) of the unknown tensor, it is possible to know a temporal change of each pattern (a product pattern, a user pattern, and a time pattern in the example of FIG. 1).

  Furthermore, when predicting related data at an arbitrary time, after the iteration is completed by the analysis end condition determination unit 27, the following predicted value calculation unit 4 is executed.

<Predicted value calculation unit 4>
The predicted value calculation unit 4 first searches the vectors T whose elements are observation times for the times T (k) and T (k + 1) before and after the predicted time t ′ designated by the input unit 1. This can be easily performed using a binary search because the elements of T are arranged in ascending order.

  When the index k corresponding to the preceding and succeeding times is found, using this, each element (i, j) of the predicted value X ′ at the time t ′ is calculated according to the following equation (14), and the predicted value X ′ is configured. To do.

Here, the predicted value X ′ is an n ₁ × n ₂ matrix, and its elements (i, j) are the elements (i, j, k) and (i, j) of the main variable X obtained as a result of optimization. , k + 1) is internally divided at a ratio of T (k + 1) -t ': t'-T (k).

In the case of t ′ <T (1), the element (i, j, 1) of the main variable X is set as the element (i, j) of the predicted value X ′. When t ′> T (n ₃ ), the element (i, j, n ₃ ) of the main variable X is set as the element (i, j) of the predicted value X ′.

<Output unit 5>
The output unit 5 outputs a tensor X, an error ε, a predicted value X ′, and pattern bases U ₁ , U ₂ , U ₃ that are main variables obtained as a result of optimization.

<Operation of time-related data analyzer>
Next, the operation of the time-related data analysis apparatus according to this embodiment will be described. First, when various values are input to the time-related data analysis device via the input unit 1, the time-related data analysis processing routine shown in FIG. 5 is executed in the time-related data analysis device.

First, in step 100, the data size n ₁ , n ₂ , n ₃ , N, M, time label T, observation vector y, observation index I, J, K, unobserved value input by the input unit 1 are input. Index I ^c , J ^c , K ^c , regularization constant λ, regularization balance constant γ ₁ , γ ₂ , γ ₃ , step size η, maximum number of iterations L, allowable error upper limit E, prediction time t ′ And stored in the storage unit 3.

In step 102, the main variable X, auxiliary variables Z ₁ , Z ₂ , Z ₃ , pattern bases U ₁ , U ₂ , U ₃ , and Lagrange multipliers A ₁ , A ₂ , A stored in the storage unit 3 are stored. Initialize ₃ to zero. In step 104, an initial value 1 is set to a variable l indicating the number of iterations.

In the next step 106, the main variable X stored in the storage unit 3 is updated. In step 108, auxiliary variables Z ₁ , Z ₂ , Z ₃ stored in the storage unit 3 and the pattern base U ₁ , U ₂ , U ₃ , and relationship patterns V ₁ , V ₂ , V ₃ are updated.

In step 110, a process of updating the Lagrange multipliers A ₁ , A ₂ , A ₃ stored in the storage unit 3 is performed. In step 112, an error ε of the KKT condition is calculated.

  In the next step 114, it is determined whether or not the error ε calculated in step 112 is larger than the upper limit E of the allowable error. When the error ε is equal to or less than the upper limit E of the allowable error, it is determined that the analysis end condition is satisfied, and the process proceeds to step 120. On the other hand, if the error ε is larger than the upper limit E of the allowable error, it is determined in step 116 whether or not the number of iterations 1 is less than the maximum number of iterations L. When the number of iterations l reaches the maximum number of iterations L, it is determined that the analysis end condition is satisfied, and the routine proceeds to step 120. On the other hand, if the number of iterations 1 is less than the maximum number of iterations L, the process proceeds to step 118, increments the number of iterations 1 by 1, and returns to step 106 above.

In step 120, based on the tensor X of the main variable stored in the storage unit 3, the relationship data at the predicted time t ′ is predicted. In step 122, the finally obtained main variable X, error ε, predicted value X ′, and pattern bases U ₁ , U ₂ , U ₃ are output, and the time-related data analysis processing routine is terminated.

  The processing in step 106 is realized by the processing routine shown in FIG.

First, in step 1060, based on the auxiliary variables Z ₁ , Z ₂ , Z ₃ and the Lagrange multipliers A ₁ , A ₂ , A ₃ , the predicted value ^ X of the main variable is constructed according to the above equation (4). In step 1062, the main variable X is updated according to the equation (5) based on the predicted value ^ X of the main variable configured in step 1060, and stored in the main variable 34 of the storage unit 3. The processing routine is terminated.

  The processing of step 108 is realized by the processing routine shown in FIG.

  First, in step 1080, the variable k indicating the mode is set to an initial value 0, and in step 1082, k is incremented by one.

In step 1084, based on the main variable X and the Lagrange multipliers A ₁ , A ₂ , A ₃ , the auxiliary variable predicted value ^ Z _k is constructed according to the above equation (6). In the next step 1086, the singular value decomposition is calculated according to the above equation (7) for the predicted value ^ Z _k of the auxiliary variable configured in step 1084. In step 1088, the singular value decomposition in step 1086 is performed. Based on the result, the pattern base U _k and the relation pattern V _k are updated and stored in the pattern base 38 and the relation pattern 39 of the storage unit 3.

Then, in step 1090, based on the matrix S _k obtained by the singular value decomposition in step 1086, the pattern base U _k updated in step 1088, and the relationship pattern V _k , the auxiliary is performed according to the above equation (8). The variable Z _k is updated and stored in the auxiliary variable 35 of the storage unit 3.

  In step 1092, it is determined whether or not the variable k is less than 3. If the variable k is less than 3, the process returns to step 1082 to increment k by 1. On the other hand, when the variable k reaches 3, it is determined that the processing of steps 1082 to 1090 has been executed for each mode k, and the processing routine is ended.

  The processing of step 110 is realized by the processing routine shown in FIG.

  First, in step 1100, a variable k indicating a mode is set to an initial value 0, and in step 1102, k is incremented by one.

Then, stored in step 1104, the main variable X, and an auxiliary coefficient Z _k, based on the Lagrange multipliers _{A k,} update the Lagrange multiplier _{A k} according to the above (9), the Lagrange multiplier 36 of the storage section 3 To do.

  In Step 1106, it is determined whether or not the variable k is less than 3. If the variable k is less than 3, the process returns to Step 1102 and increments k by 1. On the other hand, when the variable k reaches 3, it is determined that the processing in steps 1102 to 1104 has been executed for each mode k, and the processing routine is terminated.

  The processing of step 112 is realized by the processing routine shown in FIG.

First, in step 1120, it is determined whether the regularization constant λ is greater than zero. If the regularization constant λ is greater than 0, in step 1122, the above (11) is determined based on the main variable X, the observation vector y, the regularization constant λ, the step size η, and the Lagrange multipliers A ₁ , A ₂ , A _3. ), The error term ε ₁ is calculated, and the process proceeds to step 1126. On the other hand, if the regularization constant λ is 0, in step 1124, the error term ε ₁ is calculated according to the above equation (12) based on the Lagrange multipliers A ₁ , A ₂ , A ₃ , and the flow goes to step 1126. Transition.

In step 1126, the error term ε ₂ is calculated according to the above equation (13) based on the main variable X and the auxiliary variables Z ₁ , Z ₂ , Z ₃ . In the next step 1128, based on the error term ε ₁ calculated in step 112 or 1124 and the error term ε ₂ calculated in step 1126, the error ε is calculated according to the above equation (10), and the processing routine is performed. Exit.

<Experimental example>
Next, the results of an experiment using the time-related data analysis method proposed in this embodiment will be described. In the experiment, data from 1985 to 2005 was used every five years in the input-output table data (inverse matrix coefficient) published by the Statistics Bureau of the Ministry of Internal Affairs and Communications. In addition, since the original data has a period aggregated in 32 industries and a period aggregated in 34 industries, the data recalculated in 32 industries was used. Here, the data for 1995 is assumed as data to be predicted (that is, T (t ′) = 1995). Furthermore, 20% of randomly selected elements from the remaining 1985, 1990, 2000, and 2005 data were used as missing data. Under this setting, n ₁ = 32, n ₂ = 32, n ₃ = 4, N = 4096, and M = 3277. Further, the regularization constant λ = 0, the regularization balance constant γ ₁ = 0, γ ₂ = 0, γ ₃ = 1, the step size η = 1, the allowable error E = 0.001, and the maximum number of repetitions L = 1000. However, the algorithm achieved an error below the allowable error E in about 60 iterations before reaching the maximum number of iterations. Further, an error err1 for the unobserved element is defined as the following equation (15).

err1 = ‖true (I ^c , J ^c , K ^c ) −X (I ^c , J ^c , K ^c ) ‖ / ‖X (I ^c , J ^c , K ^c ) ‖ (15)

Where Xtrue is the true inverse matrix coefficient, X (I ^c , J ^c , K ^c ) is the unobserved elements of tensor X arranged in a (N−M) dimensional vector, and ‖ ‖ is the usual Represents the Euclidean norm. Further, the error err2 at the prediction point t ′ is defined as the following equation (16).

err2 = ‖Xtrue (:,:, t ') − X' (:) ‖ / ‖Xtrue (:,:, t ') ‖ (16)

  Here, Xtrue (:,:, t ') is the true inverse matrix coefficient in 1995, vertically connected in the order of the 1st column, 2nd column, ..., 32nd column, to a 1024-dimensional vector. X ′ (:) is the same as the vector of the estimated inverse matrix coefficient in 1995.

The experimental results are as follows. First, the rank of the obtained solution was r ₁ = 32, r ₂ = 32, r ₃ = 4. For mode 1 and mode 2, this is a reasonable result because the regularization balance constants γ ₁ = 0 and γ ₂ = 0. Although it is full rank regarding mode 3, it turns out that the obtained singular value concentrates on the 1st singular value as shown in FIG. Therefore, it is considered that the low rank assumption is applied. The errors obtained by repeating the random selection of the observation element (80%) 100 times were err1 = 0.11 for the unobserved element and err2 = 0.050 for the interpolation. FIG. 11 shows the true inverse matrix coefficients in 1995 and the estimated inverse matrix coefficients. Here, only the vertical axis indicates the name of the business type, but the horizontal axis also corresponds to the same business type as the vertical axis. This figure also shows that the estimation result is in good agreement with the true inverse matrix coefficient. 12 and 13 show the singular vectors associated with the singular value decomposition (FIG. 10 above) relating to mode 3. FIG. FIG. 12 shows the first component, and FIG. 13 shows the second component. 12 and 13 are multiplied by singular values corresponding to the first column (first singular vector) and the second column (second singular vector) of the pattern base U ₃ , respectively. Note that each column vector of U ₃ obtained as a result of optimization has four elements, corresponding to time points 1985, 1990, 2000, 2005, and 1995 data (circled) is 1990 And predicted values obtained from internal data of 2000 data. On the left side of FIGS. 12 and 13, singular vectors corresponding to the first and second singular values of Z ₃ , that is, the first and second columns of the relationship pattern V ₃ are transformed into 32 × 32 matrices, respectively. It corresponds to the inverse matrix coefficient component that changes with time according to the dynamics represented by the pattern base shown on the right side of the figure. FIG. 12 above shows that the linkages in the entire industry decreased year by year until 2000, and began to increase after 2000. In Figure 13 above, the impact on other industries from traditional industries such as agriculture, forestry and fisheries and steel decreases, and conversely the impact on other industries in non-traditional industries such as communications / broadcasting and services for business establishments. It shows that it is getting stronger.

  As described above, according to the time-related data analysis apparatus according to the embodiment of the present invention, the time-related data indicating the strength of the relationship between time-varying objects is analyzed to minimize the objective function. In addition, the strength of the relationship between unobserved objects can be determined by determining the unknown tensor X representing the time relationship data.

In addition, the inherent data can be found by decomposing the time-related relational data into a small number of patterns. In the singular value decomposition of the auxiliary variable update unit 24, the number of patterns is automatically optimized without setting the number of patterns in advance by processing only singular values having a value of γ _k / η or more and corresponding singular vectors. The correct value is obtained.

  It is also possible to decompose the time evolution of the relational data to be analyzed into a plurality of time factors. In addition, since the relational data to be analyzed is expressed as a weight of a plurality of patterns, it is possible to restore / predict the observation data on which unknowns or noises that have failed to be observed are highly accurate.

  Note that the present invention is not limited to the above-described embodiment, and various modifications and applications are possible without departing from the gist of the present invention.

  For example, the time-related data analysis apparatus described above has a computer system inside, but the “computer system” includes a homepage providing environment (or display environment) if a WWW system is used. Shall be.

  In the present specification, the embodiment has been described in which the program is installed in advance. However, the program can be provided by being stored in a computer-readable recording medium.

DESCRIPTION OF SYMBOLS 1 Input part 2 Analysis part 3 Storage part 4 Predicted value calculation part 5 Output part 21 Initial setting part 22 Variable calculation part 23 Main variable update part 24 Auxiliary variable update part 25 Lagrange multiplier update part 26 Error calculation part 27 Analysis end condition determination part

Claims (5)

Representing a time series of relational data indicating the strength of the relation between each of the plurality of first objects and each of the plurality of second objects, and analyzing the time relational data represented by a tensor, to analyze the time A time relational data analysis device for obtaining an unknown tensor X representing relational data,
Each mode k (k = 1) of an observation vector y consisting of M observation values of the strength of the relationship between the first object and the second object, a regularization constant λ that is a positive real number, and a tensor X , 2, 3) storing means for storing a regularization balance constant γ _k representing a relative regularization strength relative to
The observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) are acquired from the storage means, and the objective function expressed by the following equation (I) is minimized. Analyzing means for obtaining the tensor X to be converted;
A temporal relationship data analyzer comprising,
The storage means includes an auxiliary variable Z _k (k = 1, 2, 3), a value x of each element of the tensor X, a Lagrange multiplier A _k (k = 1, 2, 3), a first object direction, a second object direction, and each of the pattern base U _k in the time direction (k = 1, 2, 3), the first object that varies in accordance with said pattern base U _k (k = 1, 2, 3), the second object, And a relationship pattern V _k (k = 1, 2, 3) representing the strength of the relationship between each pair of time and
The analysis means includes
Based on the auxiliary coefficient Z _k and Lagrange multiplier A _k of each mode k, the predicted value of the tensor ^ X is calculated according to the following equation (II), and based on the calculated predicted value of the tensor ^ X and the observed vector y Main variable updating means for updating the value x of each element of the tensor X according to the following formula (III) and storing the value x in the storage means:
Based on the tensor X and the Lagrange multiplier A _k of each mode k, for each mode k, an auxiliary variable predicted value ^ Z _k is calculated according to the following equation (IV) , and the calculated auxiliary variable is predicted: Based on the value ^ Z _k , singular value decomposition is performed for each mode k according to the following equation (V) to obtain a diagonal matrix S _{k in} which singular values are arranged, and the pattern base U _k and the relationship Update the pattern V _k and store it in the storage means;
Based on the obtained diagonal matrix S _k , the pattern base U _k , and the relation pattern V _k , the auxiliary variable Z _k is updated according to the following equation (VI) and stored in the storage means. An auxiliary variable updating means for
Based on the tensor X and the auxiliary variable Z _k , for each mode k, the Lagrange multiplier update means for updating the Lagrange multiplier A _k and storing it in the storage means according to the following equation (VII) :
The tensor X stored in the storage means when the predetermined analysis end condition is satisfied by repeating the update by the main variable update means, the update by the auxiliary variable update means, and the update by the Lagrange multiplier update means. And an analysis end condition determination means for outputting
Time-related data analysis device.

Here, n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the time-related data, and x is an element of the tensor X. Ω (X) is an operation for extracting M elements corresponding to the observation value of the observation vector from the tensor X. mat ₁ is a mode 1 matrix that returns n ₁ × n ₂ n ₃ , mat ₂ is a mode 2 matrix that returns n ₂ × n ₃ n ₁ , and mat ₃ is n ₃ × n ₁ n Mode 3 matrixing that returns ₂ . ‖‖ _tr represents the linear sum of the singular values of the matrix.

However, ￣Ω (X) is an operation for extracting an element other than the observed value from the tensor X, and η is a constant. Representing a time series of relational data indicating the strength of the relation between each of the plurality of first objects and each of the plurality of second objects, and analyzing the time relational data represented by a tensor, to analyze the time A time relational data analysis device for obtaining an unknown tensor X representing relational data,
Each mode k (k = 1) of an observation vector y consisting of M observation values of the strength of the relationship between the first object and the second object, a regularization constant λ that is a positive real number, and a tensor X , 2, 3) storing means for storing a regularization balance constant γ _k representing a relative regularization strength relative to
The observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) are acquired from the storage means, and the objective function expressed by the following equation (I) is minimized. Analyzing means for obtaining the tensor X to be converted;
Based on the tensor X obtained by the analyzing means, the predicted value X of the relationship data at the prediction target time t ′ in the relationship between the first object i and the second object j according to the following equation (VIII): 'Relational data prediction means for calculating (i, j);
Time-related data analysis device including

Here, n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the time-related data, and x is an element of the tensor X. Ω (X) is an operation for extracting M elements corresponding to the observation value of the observation vector from the tensor X. mat ₁ is a mode 1 matrix that returns n ₁ × n ₂ n ₃ , mat ₂ is a mode 2 matrix that returns n ₂ × n ₃ n ₁ , and mat ₃ is n ₃ × n ₁ n Mode 3 matrixing that returns ₂ . ‖‖ _tr represents the linear sum of the singular values of the matrix.

However, the times when the observed values are observed before and after the prediction target time t ′ are T (k) and T (k + 1). An observation vector y consisting of M observations of the strength of the relationship between the first object and the second object, a regularization constant λ that is a positive real number, and each mode k (k = 1, 2) of the tensor X 3) In the time-related data analysis apparatus provided with storage means for storing the regularization balance constant γ _k representing the relative regularization strength relative to 3), each of the plurality of first objects and each of the plurality of second objects Time relationship data analysis method for representing an unknown tensor X representing the time relationship data by analyzing the time relationship data represented by a tensor and representing a time series of relationship data indicating the strength of the relationship between Because
The analysis unit obtains the observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) from the storage unit, and is expressed by the following equation (I): obtains the tensor X that minimizes the objective function that,
The storage means includes an auxiliary variable Z _k (k = 1, 2, 3), a value x of each element of the tensor X, a Lagrange multiplier A _k (k = 1, 2, 3), a first object direction, a second object direction, and each of the pattern base U _k in the time direction (k = 1, 2, 3), the first object that varies in accordance with said pattern base U _k (k = 1, 2, 3), the second object, And a relationship pattern V _k (k = 1, 2, 3) representing the strength of the relationship between each pair of time and
Obtaining the tensor X by the analyzing means
Based on the auxiliary coefficient Z _k and Lagrange multiplier A _k of each mode k, the main variable updating means calculates the predicted tensor value ^ X according to the following equation (II), and calculates the predicted tensor predicted value ^ X and Updating the value x of each element of the tensor X based on the observation vector y according to the following equation (III) and storing it in the storage means:
Based on the tensor X and the Lagrange multiplier A _k of each mode k, the auxiliary variable update means calculates the predicted value ^ Z _k of the auxiliary variable for each mode k according to the following equation (IV). On the basis of the predicted value ^ Z _k of the auxiliary variable , for each mode k, singular value decomposition is performed according to the following equation (V) to obtain a diagonal matrix S _{k in} which singular values are arranged, Updating the pattern base U _k and the relation pattern V _k and storing them in the storage means;
Based on the obtained diagonal matrix S _k , the pattern base U _k , and the relation pattern V _k , the auxiliary variable Z _k is updated according to the following equation (VI) and stored in the storage means. And steps to
The Lagrange multiplier updating means on the basis of the tensor X and auxiliary variables Z _k, for each mode k, the following (VII) wherein the step of storing in the storage means to update the Lagrange multiplier A _k,
When the analysis end condition determination unit repeats the update by the main variable update unit, the update by the auxiliary variable update unit, and the update by the Lagrange multiplier update unit, the storage unit when the predetermined analysis end condition is satisfied Outputting the tensor X stored in the time-related data analysis method.

Here, n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the time-related data, and x is an element of the tensor X. Ω (X) is an operation for extracting M elements corresponding to the observation value of the observation vector from the tensor X. mat ₁ is a mode 1 matrix that returns n ₁ × n ₂ n ₃ , mat ₂ is a mode 2 matrix that returns n ₂ × n ₃ n ₁ , and mat ₃ is n ₃ × n ₁ n Mode 3 matrixing that returns ₂ . ‖‖ _tr represents the linear sum of the singular values of the matrix.

However, ￣Ω (X) is an operation for extracting an element other than the observed value from the tensor X, and η is a constant. An observation vector y consisting of M observations of the strength of the relationship between the first object and the second object, a regularization constant λ that is a positive real number, and each mode k (k = 1, 2) of the tensor X 3) In the time-related data analysis apparatus provided with storage means for storing the regularization balance constant γ _k representing the relative regularization strength relative to 3), each of the plurality of first objects and each of the plurality of second objects Time relationship data analysis method for expressing an unknown tensor X representing the time relationship data by analyzing time relationship data expressed by a tensor and representing a time series of relationship data indicating the strength of the relationship between Because
The analysis unit obtains the observation vector y, the regularization constant λ, and the regularization balance constant γ _k (k = 1, 2, 3) from the storage unit, and is expressed by the following equation (I): obtains the tensor X that minimizes the objective function that,
Based on the tensor X obtained by the analyzing means by the relation data predicting means, the prediction target time t ′ in the relation between the first object i and the second object j according to the following equation (VIII): A temporal relational data analysis method further comprising calculating a predicted value X '(i, j) of the relational data.

Here, n ₁ , n ₂ , and n ₃ are natural numbers representing the size of the time-related data, and x is an element of the tensor X. Ω (X) is an operation for extracting M elements corresponding to the observation value of the observation vector from the tensor X. mat ₁ is a mode 1 matrix that returns n ₁ × n ₂ n ₃ , mat ₂ is a mode 2 matrix that returns n ₂ × n ₃ n ₁ , and mat ₃ is n ₃ × n ₁ n Mode 3 matrixing that returns ₂ . ‖‖ _tr represents the linear sum of the singular values of the matrix.

However, the times when the observed values are observed before and after the prediction target time t ′ are T (k) and T (k + 1). As each unit constituting the time relationship data analyzing device according to claim 1 or 2, the program for causing a computer to function.