User brian borchers - MathOverflow most recent 30 from http://mathoverflow.net 2013-05-23T13:24:11Z http://mathoverflow.net/feeds/user/9022 http://www.creativecommons.org/licenses/by-nc/2.5/rdf http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Answer by Brian Borchers for solve non-convex quadratic constrained quadratic programming Brian Borchers 2013-05-21T03:21:15Z 2013-05-21T03:21:15Z <p>It's helpful if you cite the paper in which you saw something that you're asking a question about- we could provide a better answer if we knew where the question came from. </p> <p>First, assume without loss of generality that $A$ and $C$ are symmetric matrices. It's easy to take these quadratic forms and write them in terms of symmetric matrices.</p> <p>I believe that this problem transformation requires that $A$ be a positive semidefinite matrix- see below. </p> <p>We'll begin with your second problem and show that it is equivalent to the original problem. We begin with</p> <p>$\min \mbox{tr}(AB)$</p> <p>subject to</p> <p>$\mbox{tr}(CB)=1$</p> <p>$\beta \geq 0$</p> <p>$\left[ \begin{array}{cc} 1 &amp; \beta^{T} \\ \beta &amp; B \end{array} \right] \succeq 0$</p> <p>By Schur's theorem, the constraint</p> <p>$\left[ \begin{array}{cc} 1 &amp; \beta^{T} \\ \beta &amp; B \end{array} \right] \succeq 0$</p> <p>is equivalent to </p> <p>$B \succeq \beta \beta^{T}$</p> <p>Note that we've implicitly restricted $B$ to being a symmetric matrix.</p> <p>Next, write $B$ as </p> <p>$B= \beta \beta^{T} + LL^{T}$</p> <p>where $L$ is the (slightly generalized) Cholesky factor of $B- \beta \beta^{T}$. If $B-\beta \beta^{T}$ is singular, then $L$ would be singular or even $0$. </p> <p>Then our problem is equivalent to </p> <p>$\min \mbox{tr}( A ( \beta \beta^{T} + LL^{T} ) )$</p> <p>subject to</p> <p>$\mbox{tr}(CB)=1$</p> <p>$\beta \geq 0$</p> <p>$B=\beta \beta^{T} + LL^{T}$</p> <p>By the property $\mbox{tr}(DEF)=\mbox{tr}(FDE)$, this problem is equivalent to </p> <p>$\min \mbox{tr}(\beta^{T}A\beta + L^{T}AL)$</p> <p>subject to</p> <p>$\mbox{tr}(CB)=1$</p> <p>$\beta \geq 0$</p> <p>$B=\beta \beta^{T} + LL^{T}$</p> <p>Note that because $A \succeq 0$, $\mbox{tr}(L^{T}AL) \geq 0$ for all $L$, and the term is minimized when $L=0$.</p> <p>Thus our problem is equivalent to:</p> <p>$\min \mbox{tr}(\beta^{T}A\beta)$</p> <p>subject to</p> <p>$\mbox{tr}(CB)=1$</p> <p>$\beta \geq 0$</p> <p>$B=\beta \beta^{T}$</p> <p>Since $B=\beta \beta^{T}$, and by the cyclic property of $\mbox{tr}()$, this is equivalent to</p> <p>$\min \mbox{tr}(\beta^{T}A\beta)$</p> <p>subject to</p> <p>$\beta^{T}C\beta=1$</p> <p>$\beta \geq 0$</p> <p>This was your original problem. </p> http://mathoverflow.net/questions/131200/what-does-vertex-solution-mean/131224#131224 Answer by Brian Borchers for What does "Vertex Solution" mean? Brian Borchers 2013-05-20T13:47:48Z 2013-05-20T13:47:48Z <p>In the context of linear programming, and assuming that you're using the simplex method to solve your LP's rather than an interior point method, it's most likely that the author means "basic feasible solution" (BFS) here. In geometrical terms, the basic feasible solutions of an LP are vertices of the polytope of feasible solutions. </p> <p>The number of vertices of a polytope defined by a finite system of linear equalitions and inequalities is finite and bounded by a function involving the number of variables and constraints. </p> <p>Since each iteration of the simplex method ends with a basic feasible solution, any optimal solution returned by the simplex method will be a BFS. An algorithm which solves a sequence of linear programming problems in which the constraints do not change (but the objective function varies) by using the simplex method will produce a sequence of basic feasible solutions from this finite set. </p> http://mathoverflow.net/questions/125403/minimization-of-a-function-when-the-feasible-set-is-an-unbounded-cone/125425#125425 Answer by Brian Borchers for minimization of a function when the feasible set is an unbounded cone Brian Borchers 2013-03-24T04:16:35Z 2013-03-24T04:16:35Z <p>No. The problem with this simple minded approach is that the sequence of constraints that you add might go on forever without adding a critical constraint. </p> <p>Consider the following example problem.</p> <p>$\min x$</p> <p>subject to</p> <p>$x \ge 0$</p> <p>$x \ge -1-1/n, \;\;\; n=1, 2, ...$</p> <p>Now, suppose that you start with the inequality $x \ge -2$, and that in each successive iteration you find that the next inequality of the form $x \ge -1-1/n$ is violated and add it to the problem. Your sequence of solutions will have $x$ converging to -1, which clearly isn't optimal. </p> <p>You'll need to specify additional structure in the problem and/or use a more sophisticated algorithm. </p> http://mathoverflow.net/questions/117279/projection-and-positive-matrices/117281#117281 Answer by Brian Borchers for Projection and Positive matrices Brian Borchers 2012-12-26T21:06:03Z 2012-12-26T21:06:03Z <p>What do you mean by "positive linear combination"? </p> <p>Try G=I. If you select any column of G, you'll find that its projection onto the space spanned by the other columns of G is the 0 vector. </p> http://mathoverflow.net/questions/117216/relating-the-angle-between-two-vectors-to-max-and-min-eigenvalues/117228#117228 Answer by Brian Borchers for Relating the angle between two vectors to max and min eigenvalues Brian Borchers 2012-12-26T03:32:37Z 2012-12-26T10:54:10Z <p>Here's a restatement of the problem for those who don't want to find the paper referenced in the question.</p> <p>We're interested in a symmetric and positive definitive matrix $\Sigma$ that has been orthogonally diagonalized as</p> <p>$\Sigma=Q^{T}\mbox{diag}(\theta_{1}^{2},...,\theta_{n}^{2})Q$</p> <p>where $\theta_{1}^{2} \geq \theta_{2}^{2} \geq ... \geq \theta_{n}^{2} > 0$, and we let $\theta_{\max}^{2}=\theta_{1}^{2}$ and $\theta_{\min}^{2}=\theta_{n}^{2}$. </p> <p>Also, let </p> <p>$\Theta^{2}=\mbox{diag}(\theta_{1}^{2},...,\theta_{n}^{2})$.</p> <p>We want to show for any nonzero vector $\alpha$, </p> <p>$\frac{\alpha^{T}\Sigma^{-1}\alpha}{\sqrt{\alpha^{T}\alpha}\sqrt{\alpha^{T}\Sigma^{-2}\alpha}} \geq \frac{\theta_{\max}\theta_{\min}}{(\theta_{\max}^{2}+\theta_{\min}^{2})/2}$</p> <p>Note that you can assume without loss of generality that $\alpha$ is of length 1. (just scale the length of alpha out of everything on the left hand side of the inequaility.)</p> <p>Also, by using the substitution $x=Q\alpha$, you can reduce this to a problem about the digonal matrix $\Theta$, and then reduce the matrix-vector products to sums. We then want to show for all vectors $x$ of length 1, </p> <p>$\frac{\sum_{i=1}^{n} x_{i}^{2}\theta_{i}^{-2}} {\sqrt{\sum_{i=1}^{n}x_{i}^{2}\theta_{i}^{-4}}} \geq \frac{\theta_{\max}\theta_{\min}}{(\theta_{\max}^{2}+\theta_{\min}^{2})/2}$</p> <p>Unfortunately, I don't see any easy way to proceed from here. </p> http://mathoverflow.net/questions/115706/how-to-solve-a-system-of-linear-equations-without-storing-the-matrix/115716#115716 Answer by Brian Borchers for How to solve a system of linear equations without storing the matrix? Brian Borchers 2012-12-07T14:21:30Z 2012-12-07T14:21:30Z <p>Being able to get elements of the matrix isn't very useful (particularly if you don't know where the nonzero elements of the matrix are without checking.) </p> <p>Iterative methods can be useful if you have the ability compute matrix vector products $Mx$. You haven't said whether this is possible. </p> <p>It seems quite likely that there's some special structure to your particular problem that would make it possible to simplify this computation. You haven't told us anything about where the system of equations comes from- perhaps if you explained this in some detail we could suggest ways to proceed. </p> http://mathoverflow.net/questions/110093/what-software-one-needs-to-solve-a-big-linear-system-on-a-small-computer/110094#110094 Answer by Brian Borchers for What software one needs to solve a big linear system on a small computer? Brian Borchers 2012-10-19T13:40:02Z 2012-10-19T13:40:02Z <p>Linear systems of 10,000 equations in 10,000 unknowns can easily be solved in a few seconds using double precision floating point arithmetic on typical consumer grade PC's and even laptop computers. </p> <p>If you're writing your main program in a compiled language like C or Fortran, you consider using the LAPACK and BLAS libraries to do this kind of work. </p> <p>If you're looking for a higher level language to do this, then consider MATLAB (or its open source work-alike, Octave.) Actually, Maple can do this kind of computation too, but it can be quite hard to force Maple into using floating point arithmetic rather than its normal symbolic computation mode. </p> http://mathoverflow.net/questions/106463/nonstandard-hessian-approximations-in-gauss-newton/106477#106477 Answer by Brian Borchers for Nonstandard Hessian approximations in Gauss-Newton Brian Borchers 2012-09-06T05:04:04Z 2012-09-06T05:04:04Z <p>You could certainly try this, but you'd have to do a lot of careful analysis to derive any convergence results. Among other things to consider:</p> <ol> <li><p>Your objective E(x) is more likely to have local minima due to the nonconvexity introduced by the logarithms. This could also happen with the sum of squares objective, but in practice it's uncommon for reasonably well behaved f(x). </p></li> <li><p>Your approximate Hessian will typically be fully dense, and depending on the size of $x$, this might make the solution of the equations impractical. </p></li> </ol> <p>You haven't said anything about how you're computing $J$ or how you might be able to compute the second derivatives of $f$. </p> <p>You haven't said anything about how large your vector of parameters $x$ is. If it is large, then rather than attempting this, I'd suggest using a limited memory BFGS method to avoid storing the dense $H$. </p> <p>If it's small, then depending on the difficulty of computing the second derivatives of $f$, I'd probably use the full Hessian and implement Newton's method rather than doing all of the work to show that this Gauss-Newton like method had good theoretical properties. If the second derivatives are hard to compute, then I'd just use a conventional BFGS quasi-Newton method. </p> http://mathoverflow.net/questions/106277/levenberg-marquadt-near-the-minima-for-non-zero-residual-problems/106288#106288 Answer by Brian Borchers for Levenberg-Marquadt near the minima for non-zero-residual problems Brian Borchers 2012-09-04T00:33:19Z 2012-09-04T00:33:19Z <p>The technique of multiplying the diagonal by $(1+\lambda)$ fails when you've got a 0 on the diagonal of $J^{T}J$. Having a zero on the diagonal of $J^{T}J$ can happen when you're far away from the optimal solution (in which case it tells you nothing about what might be true at optimality), or it can happen at an optimal solution (where the singularity of $J^{T}J$ can be a suggestion that the optimal solution is not unique.)</p> <p>As Thomas stated in his answer, a common approach is to use the QR factorization to solve the least squares problem at each iteration and zero out small diagonal elements of $R$. </p> <p>If you're stuck using the Cholesky factorization (e.g. if your problem is large and sparse and the Cholesky factorization can reasonably be computed but the QR factorization can't be computed), then you can avoid this problem with zeros on the diagonal by adding a small multiple of the identity to $J^{T}J$ rather than by multiplying the diagonal elements by $(1+\lambda)$. </p> <p>You should also look carefully at the scaling of your problem- bad scaling frequently leads to poor convergence in practice. </p> http://mathoverflow.net/questions/104500/cascading-minimization-problems/104502#104502 Answer by Brian Borchers for Cascading minimization problems Brian Borchers 2012-08-11T17:35:29Z 2012-08-11T17:35:29Z <p>You can solve the first LP, obtain the optimal value, and then add a constraint that the first objective has that optimal value. Then you can add your second set of constraints to the original set of constraints and solve that LP. e.g. if your original LP is </p> <p>$\max c^{T}x$</p> <p>subject to </p> <p>$Ax=b$</p> <p>$x \geq 0$</p> <p>Suppose the optimal objective value is $p^{*}$. Now suppose you want to solve a second problem</p> <p>$\max d^{T}x$</p> <p>subject to</p> <p>$Fx=g$</p> <p>$x \in U$</p> <p>This can be formulated as </p> <p>$\max d^{T}x$</p> <p>subject to </p> <p>$Ax=b$</p> <p>$c^{T}x=p^{*}$</p> <p>$Fx=g$</p> <p>$x \ge 0$</p> <p>Most simplex based LP codes can efficiently reoptimize after adding the additional constraints. </p> http://mathoverflow.net/questions/104484/linear-programming-with-or-restrictions/104487#104487 Answer by Brian Borchers for linear programming with OR restrictions Brian Borchers 2012-08-11T14:26:01Z 2012-08-11T14:26:01Z <p>Your additional constraints make the feasible region for your problem nonconvex, and thus it cannot be represented as a linear programming problem.</p> <p>If you can obtain an upper bound $M_{i}$ on the maximum value of $x_{i}$, then there is a standard approach to these problems in which the problem is formulated as a 0-1 mixed integer linear programming problem. </p> <p>First, we introduce 0-1 variables $y_{i}$, and add the constraints </p> <p>$x_{i} \geq k_{i}y_{i}$</p> <p>If $y_{i}$ takes on a 1 value, then $x_{i}$ is forced to be greater than or equal to $k_{i}$. If $y_{i}$ is 0, then this constraint is vacuous. </p> <p>Next, we add the constraints </p> <p>$x_{i} \leq M_{i}y_{i}$</p> <p>If $y_{i}$ is 0, this forces $x_{i}=0$. If $y_{i}=1$, then this constraint does nothing. </p> <p>This combination of constraints means that either $x_{i}=0$ and $y_{i}=0$ or $x_{i}\geq k_{i}$ and $y_{i}=1$. </p> <p>There are many approaches to solving the resulting 0-1 mixed integer linear programming including branch and bound methods and cutting plane algorithms. In practice, the most powerful methods (implemented in closed source commercial codes such as IBM's CPLEX as well as a number of open source noncommercial software packages) combine these two general approaches into a "branch and cut" approach. </p> <p>It is also possible to directly implement these "semi-continuous variables" within a branch and cut algorithm, and this approach does not require an upper bound $M_{i}$. This feature is available for example in IBM's CPLEX package. </p> http://mathoverflow.net/questions/103018/characterising-semi-definite-positiveness-on-vectors-with-non-negative-entries/103024#103024 Answer by Brian Borchers for Characterising semi-definite positiveness on vectors with non-negative entries Brian Borchers 2012-07-24T17:10:28Z 2012-07-24T17:10:28Z <p>Your cone $C$ is the cone of copositive matrices. The dual of C is the cone of compeltely positive matrices. See e.g.</p> <p><a href="http://mathworld.wolfram.com/CopositiveMatrix.html" rel="nofollow">http://mathworld.wolfram.com/CopositiveMatrix.html</a></p> http://mathoverflow.net/questions/102389/ease-of-calculation-of-norm/102405#102405 Answer by Brian Borchers for Ease of calculation of norm Brian Borchers 2012-07-16T23:45:09Z 2012-07-16T23:45:09Z <p>I'm assuming that by $\| A^{1/2}b-z \|$, you're referring to the 2-norm and that by $A^{1/2}$, you're referring to the unique symmetric matrix square root of $A$. </p> <p>If you can precompute $A^{1/2}z$, then you can quickly compute $\|A^{1/2}b-z\|$.</p> <p>$\| A^{1/2}b-z \|_{2}^{2}=(A^{1/2}b-z)^{T}(A^{1/2}b-z)$</p> <p>$\| A^{1/2}b-z \|_{2}^{2}=b^{T}Ab-2b^{T}A^{1/2}z + z^{T}z$</p> <p>$\| A^{1/2}b-z \|_{2}=\sqrt{b^{T}Ab-2b^{T}A^{1/2}z + z^{T}z}$</p> <p>Note that in many cases, the Cholesky factorization of $A$ can be used in place of the symmetric matrix square root. </p> http://mathoverflow.net/questions/101615/sdp-algorithms-maximally-complementary-solutions/101783#101783 Answer by Brian Borchers for SDP Algorithms/ maximally complementary solutions Brian Borchers 2012-07-09T16:03:35Z 2012-07-09T16:03:35Z <p>The standard approach is to embed the original SDP in a self dual formulation that has strictly feasible primal/dual solutions, solve the self dual formulation and then reach conclusions about the original problem from the solution of the self-dual problem.<br> See for example: </p> <p><a href="http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.37.679&amp;rep=repl&amp;type=pdf%20%22Infeasible-start%20semidefinite%20programming%20algorithms%20via%20self-dual%20embeddings%22" rel="nofollow">http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.37.679&amp;rep=rep1&amp;type=pdf Infeasible-start semidefinite programming algorithms via self-dual embeddings</a></p> <p>Unfortunately, I believe that this will only tell you when the problem has no strictly complementary primal-dual solutions and not necessarily give you a maximally complementary solution. </p> http://mathoverflow.net/questions/101130/on-random-vectors-and-eigenvectors-of-symmetric-matrices/101198#101198 Answer by Brian Borchers for On Random Vectors and Eigenvectors of Symmetric Matrices Brian Borchers 2012-07-03T01:58:57Z 2012-07-03T01:58:57Z <p>In this case, it's easy enough to compute the probability exactly in terms of the incomplete beta function.</p> <p>Let $v$ be the fixed unit vector, and $r$ be the random unit vector, uniformly distributed over the $n$-dimensional hypersphere. </p> <p>$P(| r^{T} v | &lt; \alpha) = 1-P(| r^{T}v | \geq \alpha)$</p> <p>$P(| r^{T} v | &lt; \alpha) = 1-2P(r^{T}v \geq \alpha)$</p> <p>Let $\theta$ be the angle between $r$ and $v$. Then $\cos \theta=r^{T}v$. Now, $r^{T}v \geq \alpha$ only if $r$ is in the spherical cap of the unit hypersphere centered around $v$, containing vectors within an angle $\cos^{-1}\alpha$ from $v$. The height of this spherical cap is $h=1-\alpha$. </p> <p>The probability that $r^{T}v \geq \alpha$ is then given by the ratio of the surface area of this circular cap to the surface area of the hypersphere. Using standard formulas for the surface area of the hypersphere and of the spherical cap, we get that </p> <p>$P(| r^{T}v | &lt; \alpha)=1-I_{2h-h^{2}}(\frac{n-1}{2},\frac{1}{2})$</p> <p>A quick test with a Monte Carlo simulation in MATLAB verifies this formula. </p> http://mathoverflow.net/questions/100652/efficient-algorithm-for-projection-onto-a-convex-set/100703#100703 Answer by Brian Borchers for Efficient Algorithm For Projection Onto A Convex Set Brian Borchers 2012-06-26T17:25:57Z 2012-06-26T17:25:57Z <p>You haven't told us anything about the size of your problem instances. How many terms are there in the sum of norms? What is $n$? </p> <p>Your problem is an example of a "sum of norms" optimization problem. Searching with Google Scholar will lead you to published research on this class of problems. </p> <p>CVX is using a standard approach for solving this problem by reformulating it as a second order cone programming (SOCP) problem and then using a primal-dual interior point method (SeDuMi's or SDPT3) to solve the resulting SOCP. For small problem instances, this should be a very robust and reasonably fast approach to solving the problem, but there are faster primal-dual codes for SOCP available (both CPLEX and MOSEK can be used to solve SOCP's.) You could ask CVX to extract the SOCP problem and export it in SeDuMi format and then try to use another solver on the SOCP. </p> <p>You might also look at first order methods to solve the SOCP- For example, I believe that TFOCS could be used to solve the problem. </p> http://mathoverflow.net/questions/99160/moore-penrose-bound-question/99177#99177 Answer by Brian Borchers for Moore-Penrose bound question Brian Borchers 2012-06-09T13:26:50Z 2012-06-09T23:32:48Z <p>If $Ax=b$ has a unique solution $x^{*}$, then $x_{m}=x^{*}$. </p> <p>If $Ax=b$ has infinitely many solutions, then $x_{m}$ will be one of these solutions. In particular, it will be the solution with the smallest two-norm. There will be other solutions $x$ such $\| x-x_{m} \|$ is arbitrarily large. </p> <p>If $Ax=b$ has no solutions then it isn't clear at all what you mean by referring to a specific $x$. </p> <p>Are you asking about what happens in the presence of noise in $A$ or the right hand side vector $b$? In these cases the condition number of $A$ can be used to bound the effect of the noise on the pseudoinverse solution. Unfortunately, if $\mbox{cond}(A)=\infty$, there is no bound. </p> http://mathoverflow.net/questions/98698/convolutive-noise-removal/98700#98700 Answer by Brian Borchers for Convolutive noise removal Brian Borchers 2012-06-03T04:16:31Z 2012-06-03T04:16:31Z <p>A lot depends on $\hat{\eta}(\xi)$. When you convolve this with $\hat{u}$ and $l(\xi)$, you will lose the sparsity if $\hat{\eta}(\xi)$ has broad support. </p> <p>Just how much do you know about the spectrum of $\eta(t)$? If you know it well enough, you can deconvolve it before trying your basis pursuit approach. If you don't know it very well, then you're in deep trouble. </p> http://mathoverflow.net/questions/98384/my-overdetermined-linear-system-gives-both-bad-and-good-estimates-why/98408#98408 Answer by Brian Borchers for My overdetermined linear system gives both bad and good estimates. Why ? Brian Borchers 2012-05-30T19:16:49Z 2012-05-30T19:16:49Z <p>In using least squares, you normally want the residuals from the various equations to be indpendent of each other. If your $q_{i}$ vectors are imprecise measurements, then this will introduce correlation between the residuals in the pair of equations invovling $q_{i}$. </p> <p>When you say that you get a "good fit" using the homoegneous system, how well does the solution satisfy the original nonhomogeneous system of equations? </p> http://mathoverflow.net/questions/98196/multiple-linear-regression-estimation-without-full-recalc/98204#98204 Answer by Brian Borchers for Multiple Linear Regression Estimation without full recalc Brian Borchers 2012-05-28T16:39:55Z 2012-05-28T16:39:55Z <p>There's a very large literature on updating solutions to least squares problems as new data are added. The naive formula $\hat{\beta}=(X^{T}X)^{-1}X^{T}y$ can be problematic in practice because of numerical issues with ill-conditioning. The QR factorization is typically a much better choice in practice. There are lots of papers on updating the QR factorization as data are added. See for example:</p> <p><a href="http://eprints.ma.man.ac.uk/1192/01/covered/MIMS_ep2008_111.pdf" rel="nofollow">http://eprints.ma.man.ac.uk/1192/01/covered/MIMS_ep2008_111.pdf</a></p> http://mathoverflow.net/questions/98162/a-question-about-normalization-in-the-fourier-transform/98163#98163 Answer by Brian Borchers for A question about normalization in the Fourier transform Brian Borchers 2012-05-28T03:51:15Z 2012-05-28T03:51:15Z <p>There are many variations of the definition of the Fourier transform and its inverse that are all essentially equivalent. For example, the factor of $2\pi$ often appears as $1/(2\pi)$ in front of the integral in one of the transform pair, or as $1/\sqrt{2\pi}$ in front of both integrals in the transform pair. </p> <p>Researchers in different fields have often adopted a particular convention. If you want to publish in a journal where such a convention is well established, then you will probably find that it is easier to get your paper accepted for publication if you follow the convention used by people working in that field. </p> http://mathoverflow.net/questions/97585/nonlinearly-constrained-optimization-quadratic/97590#97590 Answer by Brian Borchers for Nonlinearly constrained optimization (quadratic) Brian Borchers 2012-05-21T18:47:37Z 2012-05-21T18:47:37Z <p>The real issue here is the constraint</p> <p>$\sum_{i} x_{i}1_{x_{i}>a} &lt; b$</p> <p>whose left hand side has horrible discontinuities. </p> <p>Rather than using a solver designed for problems with continuous variables, you should formulate this as 0-1 mixed integer nonlinear programming problem, with binary decision variables $z_{i}$ that are 0 when $x_{i} \leq a$ and 1 when $x_{i}>a$. You can then use branch and bound to solve the problem. </p> http://mathoverflow.net/questions/97385/continuity-of-lexicographic-minimum-solution-of-a-parametrized-lp-problem/97436#97436 Answer by Brian Borchers for Continuity of Lexicographic Minimum Solution of a parametrized LP problem Brian Borchers 2012-05-19T21:28:13Z 2012-05-19T21:28:13Z <p>You've chosen a somewhat limited form of "parametric LP" here. It's more common to allow $F$ to vary as $F=F_{0}+tF_{1}$. For that more general variety of parametric LP, the conjecture is false- it's easy to construct a counterexample. For example, consider</p> <p>$\min tx_{1}+x_{2}$</p> <p>$x_{1}+x_{2}=1$</p> <p>$x \geq 0$</p> <p>What happens as t varies from $t&lt;1$ to $t=1$ to $t>1$? </p> <p>You could also consider parametric changes in $A$- the conjecture is also false for that case. </p> <p>Returning to the problem as stated:</p> <p>It's easy to construct examples where the LP becomes infeasible at some values of $t$. You'll need to exclude such values of $t$.</p> <p>Since this posting looks like it could well be a homework exercise, I'll simply outline an approach that you might use to prove this. First note that the basic feasible solutions (BFS's) for your parametric system of inequalities are (except at values of $t$ where a BFS vanishes or first appears) continuous functions of $t$. Also note that for each value of $t$, the set of optimal solutions is the convex hull of these BFS's. Finally, show that $x_{L}(t)$ is always at a BFS. Now, go back and think about what happens at those values of $t$ where BFS's vanish or appear. </p> http://mathoverflow.net/questions/96000/solving-for-an-operator-by-minimization/96012#96012 Answer by Brian Borchers for Solving for an operator by minimization Brian Borchers 2012-05-04T19:06:25Z 2012-05-04T19:06:25Z <p>You haven't said so, but I'm assuming that $\psi$ and $\phi$ are vectors. These could more generally be functions in some function space, and you would typically discretize those functions to work with vectors in your numerical computations. </p> <p>I'm also assuming that you have one or more $\psi$, $\phi$ examples containing measured values or values that have been computed by some numerical model. Call these examples $\psi_{1}$, $\phi_{1}$, $\psi_{2}$, $\phi_{2}$, $\ldots$, $\psi_{n}$,$\phi_{n}$. </p> <p>There may well be no exact solution to $O(x,y)\psi=\phi$ that works for all of these examples, so a typical approach would be to turn this into a nonlinear least squares problem of the form</p> <p>$\min_{x,y} \sum_{k=1}^{n} \| O(x,y)\psi_{k}-\phi_{k} \|_{2}^{2}$</p> <p>Your data and parameters are complex valued, but these are most easily dealt with by splitting the real and imaginary parts to produce a nonlinear least squares problems involving only real parameters and measured values. </p> <p>Once you've reduced the problem to a conventional real nonlinear least squares problem, you can use standard methods to solve it. The Levenberg-Marquadt method is most commonly used in practice. In the LM method, you can use finite difference approximations to get the required derivatives- most software for doing this has the ability to do this finite differencing for you. </p> http://mathoverflow.net/questions/95097/explicit-formula-for-cholesky-factorization-in-a-special-case/95100#95100 Answer by Brian Borchers for Explicit formula for Cholesky factorization in a special case Brian Borchers 2012-04-25T01:23:59Z 2012-04-25T01:23:59Z <p>The matrix $\alpha J$ is a rank one matrix, so there are simple update/downdate formulas for computing the Choleksy factorization of $Q+sI-\alpha J$ if you start with the factorization of $Q+sI$. </p> <p>I'm not aware of any update formulas that get you from the Cholesky factorization of $Q$ to a Cholesky factorization of $Q+sI$. </p> http://mathoverflow.net/questions/94960/solving-a-particular-nonlinear-system-of-equalities/94989#94989 Answer by Brian Borchers for Solving a particular nonlinear system of equalities Brian Borchers 2012-04-24T00:18:17Z 2012-04-24T00:18:17Z <p>Your problem might be small enough that it is within the range of polynomial optimization techniques based on SDP relaxations of sums of squares problem. This has been implemented in software packages such as GLOPTIPOLY. See</p> <p><a href="http://homepages.laas.fr/henrion/software/gloptipoly3/" rel="nofollow">http://homepages.laas.fr/henrion/software/gloptipoly3/</a></p> http://mathoverflow.net/questions/94459/ways-to-convert-a-positive-semi-definite-psd-matrix-positive-definite-matrix/94531#94531 Answer by Brian Borchers for Ways to convert a Positive Semi-Definite (PSD) matrix -> Positive Definite matrix Brian Borchers 2012-04-19T14:00:59Z 2012-04-19T14:00:59Z <p>The choice that you make can result in a huge difference in the solution. Neither method is particularly good and both can be quite unstable. </p> <p>For example, suppose that </p> <p>$M=[1\;0\;0\; 0;\; 0 \; 1 \; 0 \; 0; \; 0 \; 0 \; 0 \; 0;\; 0 \; 0 \; 0 \; 0]$</p> <p>and </p> <p>$N=[1;\; 1; \; 1; \; 1]$.</p> <p>What do you want $L$ to be in this situation? Why? </p> <p>Now, suppose that you try using </p> <p>$\delta = 1 \times 10^{-14}$.</p> <p>Then your approximation to $M^{1/2}$ is </p> <p>$M^{1/2}=[1\;0\;0\; 0;\; 0 \; 1 \; 0 \; 0; \; 0 \; 0 \; 1 \times 10^{-7} \; 0;\; 0 \; 0 \; 0 \; 1 \times 10^{-7}]$</p> <p>and your approximation to $M^{-1/2}$ is </p> <p>$M^{-1/2}=[1\;0\;0\; 0;\; 0 \; 1 \; 0 \; 0; \; 0 \; 0 \; 1 \times 10^{7} \; 0;\; 0 \; 0 \; 0 \; 1 \times 10^{7}]$.</p> <p>In this case, the least squares solution is </p> <p>$L=[1; \; 1; \; 1 \times 10^{14} \; 1 \times 10^{14}]$.</p> <p>In the same example, if you use the pseudoinverse, then you'd get</p> <p>$L=[1; \; 1; \; 0 \; 0]$.</p> <p>Using the pseudoinverse effectively eliminated from consideration any projection of $N$ onto the null space of $M$. That might be what you want to do. </p> <p>However, if $M$ is ill-conditioned, then it becomes practically impossible to distinguish small eigenvalues of $M$ from $0$ eigenvalues. As a result the computation using the pseudoinverse also becomes extremely unstable and dependent on the tolerance used in computing the pseudoinverse. This isn't good either. </p> <p>Continuing the example, suppose that </p> <p>$M=[1\;0\;0\; 0;\; 0 \; 1 \; 0 \; 0; \; 0 \; 0 \; 1 \times 10^{-14} \; 0;\; 0 \; 0 \; 0 \; 1 \times 10^{-14}]$.</p> <p>The pseudoinverse solution is then </p> <p>$L=[1; \; 1; \; 1 \times 10^{14};\; 1 \times 10^{14}]$</p> <p>Now, suppose that due to a rounding tolerance in computing the pseudoinverse, you compute</p> <p>$\mbox{pinv}(M)=[1\;0\;0\; 0;\; 0 \; 1 \; 0 \; 0; \; 0 \; 0 \; 0 \; 0;\; 0 \; 0 \; 0 \; 0]$</p> <p>Then you'd get </p> <p>$L=[1; \; 1; \; 0 \; 0]$</p> <p>again. </p> <p>You really need to back up a bit and explain to us why you're trying to solve this particular least squares problem in a situation where $M$ might be singular. It's likely that you need to deal with this possible ill-conditioning at a higher level by using some sort of regularization approach that is based on the particular aspects of your underlying problem. </p> http://mathoverflow.net/questions/94420/sherali-adams-relaxation/94423#94423 Answer by Brian Borchers for Sherali-Adams relaxation Brian Borchers 2012-04-18T16:19:20Z 2012-04-18T16:24:46Z <p>It isn't clear from your posting whether you're trying to understand:</p> <p>Why the inequalities generated by the Sherali-Adams procedure are valid?</p> <p>or </p> <p>Why the procedure is complete in the sense that after enough iterations you arrive at the convex hull of the integer solutions of the original integer linear programming problem?</p> <p>I suppose that you might also be interested in the question of how this can be used in practice. </p> <p>You should be aware that Sherali-Adams is just one of several "lift and project" schemes based on LP or SDP relaxations that provide a ladder of relaxations of an integer programming problem from the simple LP relaxation up to a relaxation that has only integer solutions as its extreme points. </p> <p>You might find the following paper by Monique Laurent to be useful in understanding this stuff:</p> <p><a href="http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.18.2521" rel="nofollow">http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.18.2521</a></p> http://mathoverflow.net/questions/94198/sparsity-of-qr-decomposition/94212#94212 Answer by Brian Borchers for sparsity of QR decomposition Brian Borchers 2012-04-16T13:33:12Z 2012-04-17T03:16:11Z <p>Unfortunately you can't. With any orthogonal factorization (e.g. QR, LQ, or SVD) you have the problem that because some of the columns of the orthogonal matrix have to span a particular subspace, and because the remaining columns have to form an orthogonal basis for the complement to this subspace, and because these spaces can be completely arbitrary, the orthogonal matrix won't be sparse unless you happen to be very lucky. </p> <p>There are "sparse QR" methods that effectively represent Q as a product of Givens rotations rather than storing Q explicitly. The Givens rotations are transformations that are extremely sparse/structured matrices. </p> http://mathoverflow.net/questions/93869/solving-a-non-convex-quadratically-constrained-quadratic-program/93876#93876 Answer by Brian Borchers for Solving a non-convex quadratically constrained quadratic program Brian Borchers 2012-04-12T15:49:34Z 2012-04-12T16:09:38Z <p>Branch-and-Bound can be used to solve problems like this, but it's an exponential time algorithm that is not practical for large instances of the problem. What is $n$, the dimension of the $x$ vector? How big is your $m$? </p> http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Comment by Brian Borchers Brian Borchers 2013-05-23T00:25:15Z 2013-05-23T00:25:15Z The general result is that if you have a quadratically constrained quadratic programming problem with only one non-convex constraint, then you can formulate it exactly as an SDP. That's why I believe that this relaxation will be exact in this case. However, I haven't looked into this for several years, so the details escape me at the moment. If I get a chance I'll try to sit down and work through the details. http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Comment by Brian Borchers Brian Borchers 2013-05-21T14:20:39Z 2013-05-21T14:20:39Z See appendix B of &quot;Convex Optimization&quot; by Boyd and Vandenberghe. http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Comment by Brian Borchers Brian Borchers 2013-05-21T04:02:31Z 2013-05-21T04:02:31Z On second thought, although approach I sketched out is sufficient to show that the SDP is a relaxation of the original problem, it's not enough to show that the relaxation is exact (I forgot to account for the $L^{T}CL$ term in the constraint. I'm afraid that you need the machinery of the S-lemma to show that the relaxation is exact. http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Comment by Brian Borchers Brian Borchers 2013-05-21T03:35:35Z 2013-05-21T03:35:35Z As I recall, you can also get to this formulation using the S-lemma. http://mathoverflow.net/questions/131297/solve-non-convex-quadratic-constrained-quadratic-programming/131302#131302 Comment by Brian Borchers Brian Borchers 2013-05-21T03:23:07Z 2013-05-21T03:23:07Z There are several software packages for solving such an SDP, including SDPA (written in C++) and CSDP (in C) http://mathoverflow.net/questions/123581/proof-that-polynomial-evaluated-at-roots-of-unity-is-dft Comment by Brian Borchers Brian Borchers 2013-03-05T04:31:53Z 2013-03-05T04:31:53Z Your notation seems a bit confused (and perhaps suggests why you're unable to establish this result.) You haven't explained what $a_{0}$, $a_{1}$, $...$ are. http://mathoverflow.net/questions/120634/robust-optimization-in-matlab-using-fmincon Comment by Brian Borchers Brian Borchers 2013-02-03T04:19:21Z 2013-02-03T04:19:21Z The poster has also put this same question up on scicomp.stackexchange.com, where she's somewhat more likely to get a useful answer, provided that she can clarify the question. http://mathoverflow.net/questions/117216/relating-the-angle-between-two-vectors-to-max-and-min-eigenvalues/117228#117228 Comment by Brian Borchers Brian Borchers 2012-12-26T03:44:33Z 2012-12-26T03:44:33Z You can think of the $x_{i}^{2}$ as nonnegative weights that sum to one. This opens up the whole world of the generalized mean inequality. http://mathoverflow.net/questions/115706/how-to-solve-a-system-of-linear-equations-without-storing-the-matrix/115716#115716 Comment by Brian Borchers Brian Borchers 2012-12-16T14:25:24Z 2012-12-16T14:25:24Z If the matrix isn't sparse, and the cost of getting individual matrix entries is large compared to the cost of accessing an element of a matrix stored in conventional dense matrix form, then iterative methods are going to be horribly slow in practice. http://mathoverflow.net/questions/115706/how-to-solve-a-system-of-linear-equations-without-storing-the-matrix/115716#115716 Comment by Brian Borchers Brian Borchers 2012-12-16T14:22:25Z 2012-12-16T14:22:25Z Let me clarify what I meant here- &quot;being able to get an arbitrary element M(i,j) at little cost&quot; isn't very useful. If you don't know where the nonzero elements are in the matrix, then you have to check every single one to find the nonzeros. If you do happen to know where the nonzero elements are, and you can compute them quickly, then you could use this as a way to do matrix vector multiplications in an iterative method. http://mathoverflow.net/questions/109343/optimal-gear-trains Comment by Brian Borchers Brian Borchers 2012-10-11T03:42:57Z 2012-10-11T03:42:57Z See <a href="http://www.ams.org/samplings/feature-column/fcarc-stern-brocot" rel="nofollow">ams.org/samplings/feature-column/&hellip;</a> for a discussion of the Stern-Brocot tree which has been used to design optimal gear trains. Tables of optimal gear trains based on this method have been around for a long time. http://mathoverflow.net/questions/106463/nonstandard-hessian-approximations-in-gauss-newton/106477#106477 Comment by Brian Borchers Brian Borchers 2012-09-07T13:08:16Z 2012-09-07T13:08:16Z On the other hand, if this really an unconstrained problem, I'd simply use BFGS. The computations inside the BFGS method won't be time consuming for a small problem like this, so the big concern will be the cost of function and derivative computations. You'll probably want to experiment with finite difference derivatives vs. derivatives by automatic differentiation vs. analytical formulas for the derivatives to see what works best (both in terms of accuracy and speed.) You can also play with the convergence criteria- perhaps a relatively imprecise solution will be adequate for your work. http://mathoverflow.net/questions/106463/nonstandard-hessian-approximations-in-gauss-newton/106477#106477 Comment by Brian Borchers Brian Borchers 2012-09-07T04:18:26Z 2012-09-07T04:18:26Z Now I'm even more confused about your problem. You say that x lives &quot;in a non-euclidean manifold&quot;, so it appears that you don't have an unconstrained optimization problem. In that case, you'd have a constrained nonlinear optimization problem and you'd need to consider algorithms for constrained problems rather than methods for unconstrained optimization. A good (but somewhat old) textbook that deals with mathematics and many practical aspects of optimization is &quot;Practical Optimization&quot; by Gill, Murray, and Wright. You'd find lots of good advice in that book on how to approach this. http://mathoverflow.net/questions/106277/levenberg-marquadt-near-the-minima-for-non-zero-residual-problems/106281#106281 Comment by Brian Borchers Brian Borchers 2012-09-05T18:01:45Z 2012-09-05T18:01:45Z @Alex; The $J^{T}J$ approximation to the hessian of $c(x)$ works well when $f_{i}(x)-y_{i}$ is small, but if the residuals are large then the approximation can degrade. This is because the second order terms that go into the Hessian are dropped when we approximate it with $J^{T}J$, and those terms have $f_{i}(x)-y_{i}$ factors. http://mathoverflow.net/questions/106277/levenberg-marquadt-near-the-minima-for-non-zero-residual-problems/106288#106288 Comment by Brian Borchers Brian Borchers 2012-09-05T17:59:00Z 2012-09-05T17:59:00Z This is confusing- are you minimizing a sum of squares or are you using some other objective (such as the Huber measure)? If you're not minimizing a sum of squares than LM isn't appropriate in the first place. Assuming that you are doing nonlinear least squares, and assuming that you have a reasonable assessment of the uncertainty in the measurements (e.g. $y_{i}$ is known with uncertainty $\sigma_{i}$, then you should normalize by $c(x)=\sum ((f_{i}(x)-y_{i})/\sigma_{i})^{2}$. You should also scale the parameters $x$ so that they're all of similar magnitude.