:mod:`heat.linalg.svdtools` ================================ .. py:module:: heat.core.linalg.svdtools .. autoapi-nested-parse:: distributed hierarchical SVD Module Contents --------------- .. function:: hsvd_rank(A: heat.core.dndarray.DNDarray, maxrank: int, compute_sv: bool = False, maxmergedim: Optional[int] = None, safetyshift: int = 5, silent: bool = True) -> Union[Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, float], Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray], heat.core.dndarray.DNDarray] Hierarchical SVD (hSVD) with prescribed truncation rank `maxrank`. If A = U diag(sigma) V^T is the true SVD of A, this routine computes an approximation for U[:,:maxrank] (and sigma[:maxrank], V[:,:maxrank]). The accuracy of this approximation depends on the structure of A ("low-rank" is best) and appropriate choice of parameters. One can expect a similar outcome from this routine as for sci-kit learn's TruncatedSVD (with `algorithm='randomized'`) although a different, determinstic algorithm is applied here. Hereby, the parameters `n_components` and `n_oversamples` (sci-kit learn) roughly correspond to `maxrank` and `safetyshift` (see below). :param A: 2D-array (float32/64) of which the hSVD has to be computed. :type A: DNDarray :param maxrank: truncation rank. (This parameter corresponds to `n_components` in sci-kit learn's TruncatedSVD.) :type maxrank: int :param compute_sv: compute_sv=True implies that also Sigma and V.T are computed and returned. The default is False. :type compute_sv: bool, optional :param maxmergedim: maximal size of the concatenation matrices during the merging procedure. The default is None and results in an appropriate choice depending on the size of the local slices of A and maxrank. Too small choices for this parameter will result in failure if the maximal size of the concatenation matrices does not allow to merge at least two matrices. Too large choices for this parameter can cause memory errors if the resulting merging problem becomes too large. :type maxmergedim: int, optional :param safetyshift: Increases the actual truncation rank within the computations by a safety shift. The default is 5. (There is some similarity to `n_oversamples` in sci-kit learn's TruncatedSVD.) :type safetyshift: int, optional :param silent: silent=False implies that some information on the computations are printed. The default is True. :type silent: bool, optional :returns: if compute_sv=True: U, Sigma, V.T, a-posteriori error estimate for the reconstruction error ||A-U Sigma V^T ||_F / ||A||_F (computed according to [2] along the "true" merging tree). if compute_sv=False: U, a-posteriori error estimate :rtype: (Union[ Tuple[DNDarray, DNDarray, DNDarray, float], Tuple[DNDarray, DNDarray, DNDarray], DNDarray]) .. rubric:: Notes The size of the process local SVDs to be computed during merging is proportional to the non-split size of the input A and (maxrank + safetyshift). Therefore, conservative choice of maxrank and safetyshift is advised to avoid memory issues. Note that, as sci-kit learn's randomized SVD, this routine is different from `numpy.linalg.svd` because not all singular values and vectors are computed and even those computed may be inaccurate if the input matrix exhibts a unfavorable structure. Please note that "rank" in the context of SVD always refers to the number of singular values/vectors to compute (i.e., "rank" refers to the mathematical rank of a matrix). This is completely different from the notion of "(MPI-)rank", i.e., the ID given to a process, in a parallel MPI-application. .. seealso:: :func:`hsvd`, :func:`hsvd_rtol` .. rubric:: References [1] Iwen, Ong. A distributed and incremental SVD algorithm for agglomerative data analysis on large networks. SIAM J. Matrix Anal. Appl., 37(4), 2016. [2] Himpe, Leibner, Rave. Hierarchical approximate proper orthogonal decomposition. SIAM J. Sci. Comput., 40 (5), 2018. .. function:: hsvd_rtol(A: heat.core.dndarray.DNDarray, rtol: float, compute_sv: bool = False, maxrank: Optional[int] = None, maxmergedim: Optional[int] = None, safetyshift: int = 5, no_of_merges: Optional[int] = None, silent: bool = True) -> Union[Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, float], Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray], heat.core.dndarray.DNDarray] Hierchical SVD (hSVD) with prescribed upper bound on the relative reconstruction error. If A = U diag(sigma) V^T is the true SVD of A, this routine computes an approximation for U[:,:r] (and sigma[:r], V[:,:r]) such that the rel. reconstruction error ||A-U[:,:r] diag(sigma[:r]) V[:,:r]^T ||_F / ||A||_F does not exceed rtol. The accuracy of this approximation depends on the structure of A ("low-rank" is best) and appropriate choice of parameters. This routine is similar to `hsvd_rank` with the difference that truncation is not performed after a fixed number (namly `maxrank` many) singular values but after such a number of singular values that suffice to capture a prescribed fraction of the amount of information contained in the input data (`rtol`). :param A: 2D-array (float32/64) of which the hSVD has to be computed. :type A: DNDarray :param rtol: desired upper bound on the relative reconstruction error ||A-U Sigma V^T ||_F / ||A||_F. This upper bound is processed into 'local' tolerances during the actual computations assuming the worst case scenario of a binary "merging tree"; therefore, the a-posteriori error for the relative error using the true "merging tree" (see output) may be significantly smaller than rtol. Prescription of maxrank or maxmergedim (disabled in default) can result in loss of desired precision, but can help to avoid memory issues. :type rtol: float :param compute_sv: compute_sv=True implies that also Sigma and V.T are computed and returned. The default is False. :type compute_sv: bool, optional :param no_of_merges: Maximum number of processes to be merged at each step. If no further arguments are provided (see below), this completely determines the "merging tree" and may cause memory issues. The default is None and results in a binary merging tree. Note that no_of_merges dominates maxrank and maxmergedim in the sense that at most no_of_merges processes are merged even if maxrank and maxmergedim would allow merging more processes. :type no_of_merges: int, optional :param maxrank: maximal truncation rank. The default is None. Setting at least one of maxrank and maxmergedim is recommended to avoid memory issues, but can result in loss of desired precision. Setting only maxrank (and not maxmergedim) results in an appropriate default choice for maxmergedim depending on the size of the local slices of A and the value of maxrank. :type maxrank: int, optional :param maxmergedim: maximal size of the concatenation matrices during the merging procedure. The default is None and results in an appropriate choice depending on the size of the local slices of A and maxrank. The default is None. Too small choices for this parameter will result in failure if the maximal size of the concatenation matrices does not allow to merge at least two matrices. Too large choices for this parameter can cause memory errors if the resulting merging problem becomes too large. Setting at least one of maxrank and maxmergedim is recommended to avoid memory issues, but can result in loss of desired precision. Setting only maxmergedim (and not maxrank) results in an appropriate default choice for maxrank. :type maxmergedim: int, optional :param safetyshift: Increases the actual truncation rank within the computations by a safety shift. The default is 5. :type safetyshift: int, optional :param silent: silent=False implies that some information on the computations are printed. The default is True. :type silent: bool, optional :returns: if compute_sv=True: U, Sigma, V.T, a-posteriori error estimate for the reconstruction error ||A-U Sigma V^T ||_F / ||A||_F (computed according to [2] along the "true" merging tree used in the computations). if compute_sv=False: U, a-posteriori error estimate :rtype: (Union[ Tuple[DNDarray, DNDarray, DNDarray, float], Tuple[DNDarray, DNDarray, DNDarray], DNDarray]) .. rubric:: Notes The maximum size of the process local SVDs to be computed during merging is proportional to the non-split size of the input A and (maxrank + safetyshift). Therefore, conservative choice of maxrank and safetyshift is advised to avoid memory issues. For similar reasons, prescribing only rtol and the number of processes to be merged in each step (without specifying maxrank or maxmergedim) may result in memory issues. Although prescribing maxrank is therefore strongly recommended to avoid memory issues, but may result in loss of desired precision (rtol). If this occures, a separate warning will be raised. Note that this routine is different from `numpy.linalg.svd` because not all singular values and vectors are computed and even those computed may be inaccurate if the input matrix exhibts a unfavorable structure. To avoid confusion, note that `rtol` in this routine does not have any similarity to `tol` in scikit learn's TruncatedSVD. Please note that "rank" in the context of SVD always refers to the number of singular values/vectors to compute (i.e., "rank" refers to the mathematical rank of a matrix). This is completely different from the notion of "(MPI-)rank", i.e., the ID given to a process, in a parallel MPI-application. .. seealso:: :func:`hsvd`, :func:`hsvd_rank` .. rubric:: References [1] Iwen, Ong. A distributed and incremental SVD algorithm for agglomerative data analysis on large networks. SIAM J. Matrix Anal. Appl., 37(4), 2016. [2] Himpe, Leibner, Rave. Hierarchical approximate proper orthogonal decomposition. SIAM J. Sci. Comput., 40 (5), 2018. .. function:: hsvd(A: heat.core.dndarray.DNDarray, maxrank: Optional[int] = None, maxmergedim: Optional[int] = None, rtol: Optional[float] = None, safetyshift: int = 0, no_of_merges: Optional[int] = 2, compute_sv: bool = False, silent: bool = True, warnings_off: bool = False) -> Union[Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, float], Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray], heat.core.dndarray.DNDarray] Computes an approximate truncated SVD of A utilizing a distributed hierarchical algorithm; see the references. The present function `hsvd` is a low-level routine, provides many options/parameters, but no default values, and is not recommended for usage by non-experts since conflicts arising from inappropriate parameter choice will not be caught. We strongly recommend to use the corresponding high-level functions `hsvd_rank` and `hsvd_rtol` instead. Input ------- A: DNDarray 2D-array (float32/64) of which hSVD has to be computed maxrank: int, optional truncation rank of the SVD maxmergedim: int, optional maximal size of the concatenation matrices when "merging" the local SVDs rtol: float, optional upper bound on the relative reconstruction error ||A-U Sigma V^T ||_F / ||A||_F (may deteriorate due to other parameters) safetyshift: int, optional shift that increases the actual truncation rank of the local SVDs during the computations in order to increase accuracy no_of_merges: int, optional maximum number of local SVDs to be "merged" at one step compute_sv: bool, optional determines whether to compute U, Sigma, V.T (compute_sv=True) or not (then U only) silent: bool, optional determines whether to print infos on the computations performed (silent=False) warnings_off: bool, optional switch on and off warnings that are not intended for the high-level routines based on this function .. rubric:: Notes Please note that "rank" in the context of SVD always refers to the number of singular values/vectors to compute (i.e., "rank" refers to the mathematical rank of a matrix). This is completely different from the notion of "(MPI-)rank", i.e., the ID given to a process, in a parallel MPI-application. .. rubric:: References [1] Iwen, Ong. A distributed and incremental SVD algorithm for agglomerative data analysis on large networks. SIAM J. Matrix Anal. Appl., 37(4), 2016. [2] Himpe, Leibner, Rave. Hierarchical approximate proper orthogonal decomposition. SIAM J. Sci. Comput., 40 (5), 2018. .. seealso:: :func:`hsvd_rank`, :func:`hsvd_rtol` .. function:: isvd(new_data: heat.core.dndarray.DNDarray, U_old: heat.core.dndarray.DNDarray, S_old: heat.core.dndarray.DNDarray, Vt_old: heat.core.dndarray.DNDarray, maxrank: Optional[int] = None) -> Tuple[heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray, heat.core.dndarray.DNDarray] Incremental SVD (iSVD) for the addition of new data to an existing SVD. Given the SVD of an "old" matrix, :math:`X_\textnormal{old} = `U_\textnormal{old} \cdot S_\textnormal{old} \cdot V_\textnormal{old}^T`, and additional columns :math:`N` (\"`new_data`\"), this routine computes (a possibly approximate) SVD of the extended matrix :math:`X_\textnormal{new} = [ X_\textnormal{old} | N]`. :param new_data: 2D-array (float32/64) of columns that are added to the "old" SVD. It must hold `new_data.split != 1` if `U_old.split = 0`. :type new_data: DNDarray :param U_old: U-factor of the SVD of the "old" matrix, 2D-array (float32/64). It must hold `U_old.split != 0` if `new_data.split = 1`. :type U_old: DNDarray :param S_old: Sigma-factor of the SVD of the "old" matrix, 1D-array (float32/64) :type S_old: DNDarray :param Vt_old: Transpose of V-factor of the SVD of the "old" matrix, 2D-array (float32/64) :type Vt_old: DNDarray :param maxrank: truncation rank of the SVD of the extended matrix. The default is None, i.e., no bound on the maximal rank is imposed. :type maxrank: int, optional .. rubric:: Notes Inexactness may arise due to truncation to maximal rank `maxrank` if rank of the data to be processed exceeds this rank. If you set `maxrank` to a high number (or None) in order to avoid inexactness, you may encounter memory issues. The implementation follows the approach described in Ref. [1], Sect. 2. Please note that "rank" in the context of SVD always refers to the number of singular values/vectors to compute (i.e., "rank" refers to the mathematical rank of a matrix). This is completely different from the notion of "(MPI-)rank", i.e., the ID given to a process, in a parallel MPI-application. .. rubric:: References [1] Brand, M. (2006). Fast low-rank modifications of the thin singular value decomposition. Linear algebra and its applications, 415(1), 20-30.