Difference between revisions of "SingularValueDecomp"
(13 intermediate revisions by 4 users not shown) | |||
Line 1: | Line 1: | ||
[[category:Matrix Functions]] | [[category:Matrix Functions]] | ||
− | + | == SingularValueDecomp(a, i, j, j2) == | |
− | + | [[SingularValueDecomp]] computes the singular value decomposition of a matrix. Singular value decomposition is often used with sets of equations or matrices that are singular or ill-conditioned (that is, very close to singular). It factors a matrix «a», indexed by «i» and «j», with ''[[IndexLength]](i) >= [[IndexLength]](i)'', into three matrices, ''U'', ''W'', and ''V'', such that: | |
+ | :<math>a = U W V^T</math> | ||
+ | |||
+ | where ''U'' and ''V'' are orthogonal matrices and ''W'' is a diagonal matrix. ''U'' is dimensioned by «i» and «j», ''W'' by «j» and «j2», and ''V'' by «j» and «j2». In Analytica notation: | ||
+ | |||
+ | :<code>Variable A := [[Sum]]([[Sum]](U*W, J)*[[Transpose]](V, J, J2), J2)</code> | ||
+ | |||
+ | The index «j2» must be the same size as «j» and is used to index the resulting ''W'' and ''V'' arrays. [[SingularValueDecomp]] returns an array of three elements indexed by a special system index named <code>SvdIndex</code> with each element, ''U'', ''W'', and ''V'', being a reference to the corresponding array. | ||
+ | |||
+ | To get each component in a separate global variable, use: | ||
+ | :<code>Index J2 ::= [[CopyIndex]](J)</code> | ||
+ | :<code>Variable U ::= ( , W, V ) := [[SingularValueDecomp]](A, I, J, J2)</code> | ||
+ | :<code>Variable W ::= [[ComputedBy]](U)</code> | ||
+ | :<code>Variable V ::= [[ComputedBy]](U)</code> | ||
+ | |||
+ | An less convenient but alternative approach is to use the [[Using References|# (dereference) operator]] to obtain the matrix value from each reference, as in: | ||
+ | |||
+ | :<code>Index J2 := [[CopyIndex]](J)</code> | ||
+ | :<code>Variable SvdResult := [[SingularValueDecomp]](A, I, J, J2)</code> | ||
+ | :<code>Variable U := #SvdResult[SvdIndex = 'U']</code> | ||
+ | :<code>Variable W := #SvdResult[SvdIndex = 'W']</code> | ||
+ | :<code>Variable V := #SvdResult[SvdIndex = 'V']</code> | ||
+ | |||
+ | Note that [[SingularValueDecomp]] detects whether you are capturing multiple return values or not and changes its return value to one or the other behaviors shown above. In the first case, returning each matrix as a separate return value, and in the second case returning a single return value indexed by [[SvdIndex]]. | ||
+ | |||
+ | If you need the matrices in local values, use: | ||
+ | :<code>Local (u, w, v) := [[SingularValueDecomp]](A, I, J, J2);</code> | ||
+ | |||
+ | == Matrix inverse == | ||
+ | The inverse of a square matrix A, in Analytica syntax, is | ||
+ | <code> | ||
+ | :[[Local]] Winv := [[If]] J=J2 [[Then]] 1/W [[Else]] 0; | ||
+ | :[[Transpose]]([[Sum]]([[Sum]](U*Winv, J)*[[Transpose]](V, J, J2), J2),I,J) | ||
+ | </code> | ||
+ | |||
+ | Singular value decomposition can be used for matrix inverse when the matrix A is ill-conditioned, in which case the [[Invert]] function may encounter numeric instabilities. When the matrix is ill-conditioned (the [[Determinant]] is very close to zero), then some of the elements of the diagonal of <code>W</code> will be very close to zero. To avoid the numerical instabilities, the diagonal entries corresponding to the very small <code>W</code> can be replaced with 0 in <code>wInv</code>: | ||
+ | |||
+ | <code> | ||
+ | :[[Local]] wInv := [[If]] J=J2 [[And]] [[Abs]](W)>1e-4 [[Then]] 1/W [[Else]] 0; | ||
+ | :[[Transpose]]([[Sum]]([[Sum]](U*Winv, J)*[[Transpose]](V, J, J2), J2),I,J) | ||
+ | </code> | ||
+ | |||
+ | To make this convenient to use, you can introduce a new [[User-Defined Function]] as follows: | ||
+ | |||
+ | :'''Function''' <code>MatInvert( a : [I,J] ; I,J : Index )</code> | ||
+ | :'''Definition:'''<code> | ||
+ | ::Index J2 := J; | ||
+ | ::Local (u, w, v) := SingularValueDecomp(a, I, J, J2); | ||
+ | ::Local wInv := if J=J2 And w>1e-5 then 1/w else 0; | ||
+ | ::Transpose(Sum(Sum(u*wInv, J)*Transpose(v, J, J2), J2), I, J) | ||
+ | </code> | ||
+ | |||
+ | You can then use <code>MatInvert(A,I,J)</code> in place of <code>[[Invert]](A,I,J)</code>. | ||
+ | |||
+ | == See Also == | ||
+ | * [[EigenDecomp]] | ||
+ | * [[Decompose]] | ||
+ | * [[Transpose]] | ||
+ | * [[MatrixMultiply]] | ||
+ | * [[Matrix functions]] | ||
+ | * [[:Category:Matrix Functions]] |
Latest revision as of 23:16, 23 January 2025
SingularValueDecomp(a, i, j, j2)
SingularValueDecomp computes the singular value decomposition of a matrix. Singular value decomposition is often used with sets of equations or matrices that are singular or ill-conditioned (that is, very close to singular). It factors a matrix «a», indexed by «i» and «j», with IndexLength(i) >= IndexLength(i), into three matrices, U, W, and V, such that:
- [math]\displaystyle{ a = U W V^T }[/math]
where U and V are orthogonal matrices and W is a diagonal matrix. U is dimensioned by «i» and «j», W by «j» and «j2», and V by «j» and «j2». In Analytica notation:
The index «j2» must be the same size as «j» and is used to index the resulting W and V arrays. SingularValueDecomp returns an array of three elements indexed by a special system index named SvdIndex
with each element, U, W, and V, being a reference to the corresponding array.
To get each component in a separate global variable, use:
Index J2 ::= CopyIndex(J)
Variable U ::= ( , W, V ) := SingularValueDecomp(A, I, J, J2)
Variable W ::= ComputedBy(U)
Variable V ::= ComputedBy(U)
An less convenient but alternative approach is to use the # (dereference) operator to obtain the matrix value from each reference, as in:
Index J2 := CopyIndex(J)
Variable SvdResult := SingularValueDecomp(A, I, J, J2)
Variable U := #SvdResult[SvdIndex = 'U']
Variable W := #SvdResult[SvdIndex = 'W']
Variable V := #SvdResult[SvdIndex = 'V']
Note that SingularValueDecomp detects whether you are capturing multiple return values or not and changes its return value to one or the other behaviors shown above. In the first case, returning each matrix as a separate return value, and in the second case returning a single return value indexed by SvdIndex.
If you need the matrices in local values, use:
Local (u, w, v) := SingularValueDecomp(A, I, J, J2);
Matrix inverse
The inverse of a square matrix A, in Analytica syntax, is
Singular value decomposition can be used for matrix inverse when the matrix A is ill-conditioned, in which case the Invert function may encounter numeric instabilities. When the matrix is ill-conditioned (the Determinant is very close to zero), then some of the elements of the diagonal of W
will be very close to zero. To avoid the numerical instabilities, the diagonal entries corresponding to the very small W
can be replaced with 0 in wInv
:
To make this convenient to use, you can introduce a new User-Defined Function as follows:
- Function
MatInvert( a : [I,J] ; I,J : Index )
- Definition:
- Index J2 := J;
- Local (u, w, v) := SingularValueDecomp(a, I, J, J2);
- Local wInv := if J=J2 And w>1e-5 then 1/w else 0;
- Transpose(Sum(Sum(u*wInv, J)*Transpose(v, J, J2), J2), I, J)
You can then use MatInvert(A,I,J)
in place of Invert(A,I,J)
.
Enable comment auto-refresher