We will refer to the legal set as the “feasible” set since that’s the convention. It’s very important to assume that can be encoded only using linear equations or inequalities, such a set can be equivalently called as a polyhedron (or convex polyhedron). For example, in two dimensions, regular polygons are cool whereas these guys are not cool! To summarize the simplex algorithm in one line: Dantzig observed and proved that the point we are looking for is always one of the vertices of the polyhedron, to be precise, one of the vertices optimizes or in other words, optimal solution can be assumed to be at one of the vertices. The last line is not completely true, consider the problem of or in other words, find the maximum nonnegative number. The main problem with this example is that the feasible set, in this case, the set of nonnegative numbers is unbounded. The cool fact about the simplex method is that, when it finishes (believe me that this happens in finite time),  it will either give us the point we are looking for or inform us that the problem is unbounded, that is, it can keep going in a direction and get better and better values (pretty cool or no?). Anyway, in this post we will assume that the feasible set is bounded so the answer to the LO can be found at a vertex, also known as an extreme point of the feasible set. So in principle, if you have a list of extreme points, we are done. Simplex algorithm basically computes those for us. Simplex is (and probably will be) the choice to solve LO’s even after 60 years!

To the theory folks: Simplex algorithm is not a polynomial time algorithm in the worst case, see here for an example. In fact,  many (read all) attempts to fix it were in vain. But thankfully, many smart people were able to show that LO’s can be solved in polynomial time with different algorithms, eg: Karmarkar’s algorithm, Ellipsoid algorithm. But these were so inefficient to solve practical instances, so simplex was always preferred for mysterious reasons. Recently, Daniel Spielman (and coauthors) came up with a new concept called Smoothed analysis which successfully explained the performance of the Simplex algorithm. This post is not about the Simplex algorithm but I felt like these are needed to appreciate FW method and it sets up the stage very nicely. Practically, let’s say you have purchased a commercially available LO solver eg: CPLEX, Gurobi etc., then can I use it to solve say convex optimization problems? The answer is a resounding yes provided you’re 1. smart, 2. Patient, 3. give up few things. This is similar to using a laptop (internet) to watch TV shows. Smart: know to use an internet browser, a web search engine; Patient: TV cables can transmit huge amount of data much much faster than internet as of now; 3. Give up: can’t watch all the shows (obviously, right?) but can watch the popular/important shows.

FW algorithm solves the following question: Let’s say we are given with the following two types of data viz., Here is the given data: 1. a differentiable convex function ; 2. a feasible set that is convex and compact .Now how can I pick the point in   where achieves the minimum value. Few comments are in order: 1. the class of optimizable functions is much bigger (a linear function is convex and differentiable); 2. is much bigger in some sense, bounded (so compact) polyhedrons (so convex) are cool, but we do require to be bounded and for the uninitiated, think of compact being bounded and that inequalities are always of the form or and not strict or . You can only have equations that are linear (or affine) because: is convex if and only if is affine (exercise!). Alright, here is the algorithm:

• Pick a point .
• Repeat this:
  1. Compute the gradient
  2. Compute
  3. Set

Picking  is not actually required, since step 2. will anyway make the next iterate to be in so that’s cool. That leads us to the following amazing observation: assuming that we can compute reasonably well, we have essentially converted a nonlinear (but convex) optimization into a series of optimization problems with linear function that provably converges! If is a polyhedron, then we solve a LO at each iteration. Secondly, we will see how to choose the step sizes when we see the proof. So basically the algorithm at a very high level is impressive if step 2. can be done fast, it turns out that this is true for many problems in machine learning (I’ll post some references later).

Convergence: We will just discuss the most basic convergence proof. I should mention that lot of work has gone into making the algorithm faster and the proof more tighter in general and in specific cases. Intuitively, if is linear, then we are done in one iteration by setting , surprise yay! This suggests that we need a measure of the linearity of the function, let’s call it , the curvature constant of .

Definition:

Basically it measure how about the linear approximation is for a given function on the feasible set .

Definition: The Wolfe dual is defined as where .

Weak Duality: for all .

Proof: Take any . By the minimization in the definition of , we have that,

where the second inequality is due to the convexity of .

Denote the primal-dual gap at by defined as,

So, by weak duality we have that, . With these we prove a step wise improvement in the value of the objective function.

Stepwise improvement: For any with arbitrary step size , we have that,

Proof: Let . From the definition of ,

Rewriting the above inequality and using the definition of we have that,

.

Main convergence:

Proof: Define . Use the stepwise improvement lemma and use induction. For the induction step use the step sizes . TODO: Finish the proof with the details. .

Remarks: 1. No need to use a step size strategy as in many other algorithms because of the predetermined stepsizes; 2. convergence is pretty slow in general since it is sublinear; 3. proof goes through if you choose to use line search strategies(for example, wolfe conditions) as well; 4. Need to be small or at least bounded for the convergence to work. Next time we will see some examples.

Assume be a random matrix with . The algorithm to completing the missing entries is as follows: 1. Set the entries that are missing to be ; 2. Compute the best low rank approximation  for . Then we can prove the following theorem.

Theorem: Assuming that the entries are missing at random observed with probability , we can show that where  is the best rank approximation of, observed entries and , , fill the unknown entries with .  In essence, matrix completion works if we observe  entries.

Proof: 1. , here note that the first term is always less than (or equal to) the second term since  is the best possible matrix. So we have that,

. Here  is the matrix with entries,  independent with mean . Hence by observing the fact that . This further implies that .

2. Now operator norm can be easily bounded since . So we have that , dividing  both the sides gives us the result. 

The proof as of now just consists of equations and math, I’ll try to add more explanations soon.

Matrix Bernstein Inequality (MBI): Let be  mean 0, symmetric matrices such that . Then we have that where , the constants are hidden in notation and this will be often used, basically this is analogues to the big-Oh notation.

Remark 1: This inequality is sharp for as we saw earlier.
Remark 2:

Here denotes the operator norm of the matrix, that is, for symmetric matrices, or it’s largest eigenvalue and corresponds to the “largest” eigenvector.

Proof: Note that since is a symmteric matrix, it can be written as where corresponds to the eigenvalues and eigenvectors respectively by the spectral decomposition and we assume that . Now if we have a real valued function from reals, , then for a symmetric matrix  is defined as . Alternatively, we can use matrix Taylor’s series and apply on the coefficients. We will need two important matrix inequalities which we state without proof:

1. Golden-Thompson inequality: . When are reals, we get equality (this is because reals are commutative in multiplication whereas matrices aren’t).
2. Trace inequality: Roman mentioned that this is Lieb’s inequality, but I think it’s much a simpler form than that so I’m calling it as just trace inequality, it states that: Let be a fixed symmetric matrix and consider the function , denotes the trace or the sum of the diagonal entries. Then we have that is concave on symmetric PSD matrices.

Combining both the above properties and Jensen’s inequality (), we have that, . In particular, when , we get,

which we will refer to as the Lieb’s inequality from now on. Basically this inequality allows to take the expectation inside the exponential operator, pretty neat right?

Now we ready to prove MBI. It consists essentially of two steps: 1. Markov’s inequality or union bound on the sum of random matrices; 2. Lieb’s inequality.

1. Let . Then . The inequality is due to the union bound, denotes the largest eigenvalue.Now we have the following equalities,

   I don’t know why there’s a line break, my apologies but I like this since you can see that doing and are commutative operations.

   So denoting . The last inequality is due to the fact the is nonnegative and nondecreasing.

2. Lieb’s: The trick is to deal or apply Lieb’s with one random matrix at a time and continue doing it.

It is easy to check that (or I’ll probably show this in the next post): If , then we have that, . Applying this to the result from 2 (Lieb’s) and optimizing for , we have that

.

In the next post we will apply MBI to the matrix completion problem.

The problem is very simple to state: let’s say we have a matrix of of size , this simply means we have a 2-D data array with  rows and  columns. The rows and columns are indexed over, say, users and movies respectively. The famous Netflix prize problem had a matrix like this and of course, not every user has seen every movie. So the questions is, can we fill in the entries of the matrix otherwise incomplete matrix? Indeed, they had a test dataset (or simply test set) on which your algorithms will be tested on, whoever gets the most accuracy wins the challenge, so they get lots of $$$. This problem is obviously not well posed (why?) but as with most machine learning problems, mathematical ill-posedness doesn’t really matter as far as we can do good on the test set (this is not 100% true). In fact, the name of the game is to find mathematical theory/algorithms/models that can explain the problem/solution in practice. There are actually very many cool examples of this, probably we will see about this in later posts. Anyway, to analyze matix completion problem, we will need some tools from, of course, random matrix theory (I was actually surprised that these tools were actually developed very recently, I mean within the last 10-15 years atmost, if I’m wrong, please comment). So in this post we will look at some basic concentration inequalities.

Definition of subgaussian random variable: Let be a random variable. Then, the following are equivalent (TFAE):
1. Tails:
2. Moments:
3. MGF of :
4. MGF of :

We will refer to these numbers before the statements through out this post. denotes the probability and denotes the expectation.
In these statements, s can be multiplicative factors of each other but this really doesn’t matter since we are interested in the asymptotics. Examples include: gaussian, bernoulli and nonexamples include exponential, poisson etc.. Finally we will denote smallest value of in 3 (MGF of ).

Proposition 1: If ‘s are independent subgaussian random variables with mean 0, then is also a subgaussian random variable.

Proof: We will show 4 for .

The first equality is true because of the independence of the random variables , the inequality is due to 4. Choosing to be we are done.

Hoeffding’s Inequality: Let ‘s be independent subgaussian random variables with mean 0. Then

Proof: This follows from proposition 1 and the equivalence of the definitions.

OK, now we need to know a bit about subexponential random variables.

Definition of subexponential random variable: Let be a random variable. Then, the following are equivalent (TFAE):
1′.
2′.
3′. If , then

In this case we denote  = smallest .

Bernstein’s Inequality: Let ‘s be independent subexponential random variables with mean 0. Then

Proof (sketch): This requires the use of Markov inequality on the sum and then we use 3′.

Corallary: where .

Synopsis: We have Hoeffding’s inequality for subgaussians and Bernstein’s inequality for subexponentials, so what? well, if you look at the stuff inside the exp, then you’ll notice that the subgaussians are much behaved than subexponentials since they go down (concentrate) faster or better (quadratic in subgaussians vs linear in subexponentials). In the next post we will see how this can be extended to the case where are matrices.

To answer that question, we need a little bit of machinery unfortunately, here we go. We’re going to define what it means for two 2-tuples to be linearly independent or rather linearly dependent. I’m just going to say what it means right away. Two 2-tuples are linearly dependent if we can multiply one of the two by a scalar to get the other, that is, 2(1,2) = (2,4), so (1,2) and (2,4) are linearly dependent, that’s it! A counter example will be (1,0) and (0,1) since we can’t obtain one of them from the other by multiplying any number. I cheated right there because I didn’t tell you what a scalar is. A scalar is the underlying set on which the 2-tuple is defined, that is, it can be z or q. Note that if it’s q then the set of tuple is called as a v.s. whereas if it’s z then it’s not a v.s. (it’s called as a z-module). Of course, one might immediately see in the first example that (2,4) times any integer is never equal to (1,2) and conclude that both the tuples are linearly independent, this is exactly one of the problems with z-modules. Just to finish off this, we say that a set consisting of 2-tuples is linearly independent if each tuple is linearly independent from the rest, that is, any tuple cannot be written as a scaled sum of the rest.

Observe that the 2-tuple set is much bigger than the underlying scalars be it q or z. It is much bigger than twice the size of the underlying scalars which is easy to see. We can see that by the following construction. Take a piece of paper, draw to lines that are perpendicular to each other. Mark the origin and choose a specific distance to be 1 unit. Then the horizontal and the vertical lines denote the real numbers. Now we will mark 1, 2, 3,… and so on, on both of these axes. Then z^2 consists of all the points that have integer coordinates; r^2 consists of everything on the paper. If we imagine the integers on the horizontal line to be denote z as a set of integers, we can easily see that are more than 2 such sets. To be a bit fancy, z^2 is often referred to as lattice and some might say that z^2 is embedded in r^2 (or q^2). So what we have done is enlarged our space by a huge amount just by including a coordinate. But to store any point or tuple we just need to store two numbers, that is, a tuple (3,4) can be written as 3(1,0) + 4(1,0). What’s the big deal? (1,0) and (0,1) are called as the generating set, that is, any tuple can be stored in terms of these two tuples or in other words (1,0) and (0,1) generate (or span) z^2. Note also that the same two vectors span r^2 and q^2.

We have seen two concepts viz., linear independence (or dependence) and generating set or spanning set. The key property of v.s. that z-modules don’t have is that the minimal generating set is equal to the maximal independent set. Let me explain this. Let’s say that we have a set that consists of a  2-tuple, we add a 2-tuple that is linear independent and continue to do is. In r^2, let’s say we pick (2,3) to be the first tuple, and we have added (4,5), it’s easy to see that these two 2-tuples are linearly independent because multiplying one 2-tuple by any real number will never result in the other 2-tuple. But the problem is that we can’t add any more 2-tuples because -5/2(2,3) + 3/2(4,5) = (1,0); 2(2,3) -1(4,5) = (0,1), which means that any 2-tuple can be written as a scaled sum of the 2-tuples that we already have. A fancier way of saying this is that any 2-tuple can be written as a linear combination of (2,3) and (4,5) in r^2. In other words, in r^2 these two 2-tuples are linearly independent and they span the entire space r^2. Now one can easily see that (2,3) and (4,5) are linearly independent even in the z-module as well. But the problem is that they don’t span z^2, for example, they cannot generate (1,0) anymore. So we can continue adding more elements to our set to make linearly independent sets. This means that the size of maximal independent set is not 2 anymore for z-module whereas for v.s. of 2 – dimensions, once we choose two linearly independent tuples, we can’t find any more linearly independent tuples. Though in this case we can see that we can pick (1,0) and (0,1) to generate z^2. But the point is that you can have a different linearly independent set generating z^2 and I can have a different generating set and it might not be clear how to convert one to other. (1,0) and (0,1) are also called as basis or standard basis tuples.

This phenomenon is often referred to as duality, that is, maximum of something is equal to the minimum of something else. v.s. have this duality property whereas z-modules and many others don’t have this, that is, for z-modules there is nothing called as basis tuples whereas it is easy to choose standard bases for v.s.. We can easily extend this to any finite dimensional v.s.. Now it’s easy to see why v.s. are defined over a field, stare at the examples for a while if you’re still not convinced.

Many of us take the most basic courses in college like calc, linear algebra, maybe discrete math as part of a ‘requirement’ for the degree that you intend to get in the future. ‘Requirement’ is the keyword. Imagine that you want to buy a car or bike, the first thing that one might consider doing is to learn how to drive or ride one of those. This seems to be somehow very obvious to all of us because, well because we all agree. In the same way, if one wants to get a degree in say biology, one should know how to make simple calculations like the growth of an organism, size etc.. But the problem is that the basic math courses are offered by the math department. In other words, you learn math from people who make their living out of math. So what? Now it puts an immense pressure on the teachers on how to teach these basic courses. Why, they’re after all experts duh? It turns out that in most cases, the courses are specifically designed so that if one takes these courses, they will be able to do specific aspects of computations in their fields, that is, biology or engineer or whatever. Hence as a mathematician, they are forced to dumb down a lot of things, that is, they skip lot of details, exams become hard because they don’t have enough time to cover everything in the class. The people who go to office hours and ask questions find this out very soon and manage to do well. If you do well, often you don’t hate it, it is only when you spend some time on math and you don’t well when the hate becomes real. One specific example is the Linear Programming course at UW Madison, it always excites me whenever I talk about this course. First of all, it is a requirement for any graduate student in Industrial Engineering whereas it is counted for credits for math, cs, stat and maybe many others. Well the course is an application of linear algebra basically which is a course that all the undergrads in math, stat, cs, engineering take. But the problem is that the cs, engineering don’t use linear algebra in the way that is expected of as a prerequisite for the Linear Programming course, for example, they might have never worked with words like null space, image space, span etc. in their entire life. Now this is a problem because the basic component of one of the basic algorithms called as the simplex algorithm is just entirely based on this. This is not the case with the other department students usually assuming they have taken few other math courses, because either they actually took some more linear algebra courses or because they got used to definitions. The second reason is the secret ingredient because the first reason is like any other skill, once we spend enough time on something we know about that thing in some detail (hopefully). The second reason is not a direct reason and this is how one can learn math. I can quote lot of mathematicians, for example, Neumann – ‘Young man, in mathematics you don’t understand things. You just get used to them’. This is not a really good statement in some sense and I wasn’t comfortable with this. But a more matured version of this statement (or what he really meant) can be seen from Alexander Grothendieck who said, ‘The … analogy that came to my mind is of immersing the nut in some softening liquid, and why not simply water? From time to time you rub so the liquid penetrates better, and otherwise you let time pass. The shell becomes more flexible through weeks and months—when the time is ripe, hand pressure is enough, the shell opens like a perfectly ripened avocado!‘.

Probably I should have put that quote in the beginning. Let’s assume that we would like to crack a nut. One way to do this is take a hammer and try to break it and the other way is from Grothendieck. When you take a math course to apply or use it in your field, then mathematicians are forced to teach it in the hammer way, that is, they teach math as an array of tricks that can be used to solve problems of a certain type and this array seems so artificial which either makes one wonder how they came up with it (and the one starts liking it) or one starts to hate the subject. (Again, this is not their choice or how they want to do it, it is because they are forced to this.) Also, the nut becomes something special or in other words we don’t need a hammer to break a fruit whereas when we are trying to crack a nut, it becomes extremely important that we use a hammer to open or we will never be able to see what’s inside. What we are trying to do in this blog is the other way, that is, I deliberately digress, you can think of it as surrounding the nut with water. This process is slow but eventually it will pay off in the sense that we can think of hammer as as certain type of liquid that can penetrate much more. We can also do things in a much more optimal way since we can just use the hammer at places that are hard for the liquid to penetrate. This means that the idea of liquid is much more powerful than the idea of hammer which is what we are aiming for because then we learn about the nut in a much more detail way in one shot and as usual everything is nice and beautiful.

We first choose a point in the east, preferably the places in Japan where the sun apparently rises first, let’s call the point 1. We need only one coordinate since the latitudes are not useful (perhaps there’s a fancier way to make timezones using latitudes but we will not go there). Similarly we choose the point in the west where the sun sets the last, let’s call the point 0. This means that each point (or place on the earth) corresponds to a particular number between 0 and 1. Of course lots of points share the same number. Now let’s say that we have a function or a mapping defined as f that takes any number in the interval 0, 1 and outputs a number or time in our case. What we will now do is impose some conditions on f so that it becomes mathematically nice. One nice thing to have is a continuous f and perhaps smooth too. Note that we don’t have continuity in the current world since there are uneven jumps at random points on the longitudes. So the next thing we might want is f be a nice polynomial. This is actually not required but let’s say we need this just because it helps explanation easier. That’s it!

What have we done? We replaced the arbitrary timezones with nice mathematical surrogates. We will now see the immediate consequences of this. Let me first warn you that this is not looking good. The above constructions say that whenever we move in along longitudes we need to use the polynomial to figure out the time, which means every time you travel east or west even if it’s just two blocks, time changes include your travel time plus the calculation given by the polynomial. Nightmare! Well, not so much because these numbers will be fairly insignificant. But when you actually travel 200 or 300 miles these will be significant which is nice. We will also have calculators (and apps) to do this for us so it’s all good. Another thing is that once everyone gets used to this setup, it will be easy for everyone. What does this achieve overall? First, one will be better at one of the fundamental ideas in applied math and some areas in computer science, engineering in a vague sense (not that one should be, but why say no to something free and good). Second, the world will be better off with continuous functions than otherwise (I’ll leave it for you to think about)!

A mapping also called as a function can be thought of as a machine that takes an element from a set as input and gives us (or outputs) another element from a different set (or the same set). From the definition of a mapping, it is clear that we expect only one output for every input. A relation is then just a mapping of input and out put. A relation has to satisfy three other properties which are not really important for our discussion, so we will not go over them here.

Few more words to know, aarghh! The input set is called as the ‘domain’ and the output set is called as the ‘range’. I should mention that people use different words for ‘range’, one such word is ‘codomain’ and other is ‘image’ but this is very rare. So whenever anyone uses one of these words, clarify with them what they really mean, why? Let’s look at an example to understand why it is ambiguous. Imagine that we wanna know the efficiency of a car by measuring the mileage when we drive at different speeds. As everyone knows (look it up otherwise), the efficiency varies with different speeds and there is a certain range of speeds at which the car is the most efficient, let’s call the mileage at this speed to be 20 miles per gallon and the worst performance was 10 miles per gallon. For someone who knew nothing about how a car works, the mileage of this car can be any positive number but we know that the car can never go beyond 20 miles per gallon. So the codomain or the range in this case will be all positive numbers whereas the image (or the range, this is why it’s ambiguous) will be just all numbers between 10 and 20. I prefer to use codomain for the positive numbers and range for numbers between 10 and 20 in the last example. But I might use them interchangeably, I apologize for that, but I think the meaning will be clear from the context (hopefully).

These definitions for someone looking at it for the first time may seem to be general, well in fact it is general. Mathematicians have developed tools to study about these when the domain and the codomain have some special properties i.e., (surprise!!) when the domain and the codomain are numbers. So we will for the most part assume that the domain and the codomain are numbers of some sorts. To be precise, our domain will be r^n and the range will be r^n, what the heck is this? r or R is short for ‘real’ meaning all the numbers that we know, eg: -2, square root of 2, 10, 500 etc.. There are other things as well like the set of integers denoted by Z, natural numbers N, rational numbers Q, irrational numbers I, complex numbers C etc.. Curious readers can wiki it and those who have heard of them before know that each of them are what they are. For now let’s assume that n is equal to 1 i.e. all our functions are from r^1 to r^1 (written also as r to r). What makes these sets so special? Well there are lots of reasons and in fact huge part of what mathematicians do is identifying these sets and study them. Rambling apart, these are nice sets because any two elements of these sets can be compared i.e., 5 is greater than 6 (5>6), -square root of 2 is less than +square root of 2 and so on.

One special property of R, Q (and C) is that there exists an element in R, Q (and C) respectively for any element (other than zero) such that the product of both them is 1 or the unit element. Too many quantifiers in one single line, not cool brah! What I meant was that if we pick a number from r, let’s say 14, then 1/14 times 14 is 1 and this fact is true for any nonzero number we pick, that’s it! This doesn’t work if our set is Z, because 14 multiplied with any other integer will never be equal to 1. These sets are called as fields (to be pedantic, they also should satisfy other conditions but those are obvious ones like addition, multiplication etc.). Vector spaces are sets that are defined ‘over’ these fields like Q, R, C. Intuitively, (finite dimensional) vector spaces are tuples of elements from these fields i.e. an element from our favorite 2 dimensional vector space or plane is (3,4) or in general (x,y) where we require x and y to be elements of a field. So a n – dimensional vector space over a field is just the set of all n – tuples. Yep, that’s a very big set but these are the most nicest sets that one can encounter in their mathematical life, why? We’ll see that next time and more (I promise linearity II isn’t far, hang in there!).

What is a function? For us, a function is a box or a machine which takes a number as input and  outputs ‘a’ number, pretty simple right? In our example, we say that c is a function of x, in short c = f(x). Once someone gives us the numbers a and b one can quickly compute the value c for a fixed x by using simple high school algebra techniques. Wait a minute, we can also rewrite the equation as x = c – b divided by a, tada (!) now we have a x as a function of c, x = g(c). But there’s this one small thing, we can only to do this if a is not equal to zero, ouch! We now turn to invertible functions. We say that a function y = f(x) is an invertible function if f is a function and if for every value of y there is a corresponding value of x, that is, given the output we should be able to figure out what the input should have been. When a is zero, c = b, which means no matter what x is, our output will always be b, in other words given the total cost one can never figure out how much sugar was bought.

But mathematicians are always looking for invertible functions, why? simply because we can work with y or x and convert the final answer by applying f (or the inverse of f) appropriately.  So, we want linear functions to be invertible and so we define a not be equal to zero, goodbye constant functions! It’s kind of odd that it is actually difficult to swallow this because when I was in my high school, I was told that straight lines are linear functions and constant functions are straight lines so why throw them away? Don’t be mad at me, it is what it is. To be precise, I’m trying to motivate linear (invertible) transformations which we will see next time. In fact, we actually don’t get much value allowing constant functions, we don’t have much to talk about them since they’re just constant functions (MEH). Next time we will go much deeper and start talking about matrices.