Sunday, July 15, 2012

Carrying Optimal Change in Your Wallet


The Alternative COIN-OR
How much change should I carry in my wallet? How many pennies, nickels, dimes, and quarters? If I carry too much, my pockets feel heavy, and if i carry too little, then I don't have enough to cover transactions in my office cafeteria and have to expend a dollar bill (of negligible weight) that was intended for the high-calorie vending machine, and end up carrying a lot more change after the transaction. We assume that the transaction amounts are small enough so that credit-cards are not used.

The cowboy uses a horse to cross the street whereas the O.R person uses CPLEX for the same reason. Let's analyze some simple cases here using CPLEX.

Problem 1: Determining the optimal coin distribution for a given change request.
Given that we are required to provide exact change of value v, where v = 1, ..., 100, what is the optimal distribution of coins for each scenario such that:
a) The total number of coins is minimized
b) The total weight of coins is minimized

The US Mint website provides the following coin weights in grams and its corresponding value that is captured via the following piece of java code:

static String[] name = {"penny","nickel", "dime", "quarter"};
static double[] weight = {2.500, 5.000, 2.268, 5.670};
static double[] value = {1., 5., 10., 25.};

As we can see by comparing their "bang to buck" ratios, the nickel is quite inefficient whereas the dime and quarter do pretty well. Therefore, one would expect the solution to carry more dimes and quarters, and fewer nickels and 'unscalable' pennies.  Furthermore, given that the b2b function is kinda non-convex, it's a good idea to solve these problems as discrete optimization models.


The optimal  (min count) vector of integer coin quantities (x) for a given input desired value = 'v' is determined using CPLEX-Java via the following integer knapsack model:

cplex.addMinimize(cplex.sum(x));
cplex.addEq(v, cplex.scalProd(value, x));

value count weight penny nickel dime quarter
1 1 2.5 1 0 0 0
2 2 5 2 0 0 0
3 3 7.5 3 0 0 0
4 4 10 4 0 0 0
5 1 5 0 1 0 0
6 2 7.5 1 1 0 0
7 3 10 2 1 0 0
8 4 12.5 3 1 0 0
9 5 15 4 1 0 0
10 1 2.268 0 0 1 0
11 2 4.768 1 0 1 0
12 3 7.268 2 0 1 0
13 4 9.768 3 0 1 0
14 5 12.268 4 0 1 0
15 2 7.268 0 1 1 0
16 3 9.768 1 1 1 0
17 4 12.268 2 1 1 0
18 5 14.768 3 1 1 0
19 6 17.268 4 1 1 0
20 2 4.536 0 0 2 0
21 3 7.036 1 0 2 0
22 4 9.536 2 0 2 0
23 5 12.036 3 0 2 0
24 6 14.536 4 0 2 0
25 1 5.67 0 0 0 1
26 2 8.17 1 0 0 1
27 3 10.67 2 0 0 1
28 4 13.17 3 0 0 1
29 5 15.67 4 0 0 1
30 2 10.67 0 1 0 1
31 3 13.17 1 1 0 1
32 4 15.67 2 1 0 1
33 5 18.17 3 1 0 1
34 6 20.67 4 1 0 1
35 2 7.938 0 0 1 1
36 3 10.438 1 0 1 1
37 4 12.938 2 0 1 1
38 5 15.438 3 0 1 1
39 6 17.938 4 0 1 1
40 3 12.938 0 1 1 1
41 4 15.438 1 1 1 1
42 5 17.938 2 1 1 1
43 6 20.438 3 1 1 1
44 7 22.938 4 1 1 1
45 3 10.206 0 0 2 1
46 4 12.706 1 0 2 1
47 5 15.206 2 0 2 1
48 6 17.706 3 0 2 1
49 7 20.206 4 0 2 1
50 2 11.34 0 0 0 2
51 3 13.84 1 0 0 2
52 4 16.34 2 0 0 2
53 5 18.84 3 0 0 2
54 6 21.34 4 0 0 2
55 3 16.34 0 1 0 2
56 4 18.84 1 1 0 2
57 5 21.34 2 1 0 2
58 6 23.84 3 1 0 2
59 7 26.34 4 1 0 2
60 3 13.608 0 0 1 2
61 4 16.108 1 0 1 2
62 5 18.608 2 0 1 2
63 6 21.108 3 0 1 2
64 7 23.608 4 0 1 2
65 4 18.608 0 1 1 2
66 5 21.108 1 1 1 2
67 6 23.608 2 1 1 2
68 7 26.108 3 1 1 2
69 8 28.608 4 1 1 2
70 4 15.876 0 0 2 2
71 5 18.376 1 0 2 2
72 6 20.876 2 0 2 2
73 7 23.376 3 0 2 2
74 8 25.876 4 0 2 2
75 3 17.01 0 0 0 3
76 4 19.51 1 0 0 3
77 5 22.01 2 0 0 3
78 6 24.51 3 0 0 3
79 7 27.01 4 0 0 3
80 4 22.01 0 1 0 3
81 5 24.51 1 1 0 3
82 6 27.01 2 1 0 3
83 7 29.51 3 1 0 3
84 8 32.01 4 1 0 3
85 4 19.278 0 0 1 3
86 5 21.778 1 0 1 3
87 6 24.278 2 0 1 3
88 7 26.778 3 0 1 3
89 8 29.278 4 0 1 3
90 5 24.278 0 1 1 3
91 6 26.778 1 1 1 3
92 7 29.278 2 1 1 3
93 8 31.778 3 1 1 3
94 9 34.278 4 1 1 3
95 5 21.546 0 0 2 3
96 6 24.046 1 0 2 3
97 7 26.546 2 0 2 3
98 8 29.046 3 0 2 3
99 9 31.546 4 0 2 3
100 4 22.68 0 0 0 4

b) The minimum weight solutions are generated via:

cplex.addMinimize(cplex.scalProd(weight, x));
cplex.addEq(desiredValue, cplex.scalProd(value, x));


value count weight penny nickel dime quarter
1 1 2.5 1 0 0 0
2 2 5 2 0 0 0
3 3 7.5 3 0 0 0
4 4 10 4 0 0 0
5 1 5 0 1 0 0
6 2 7.5 1 1 0 0
7 3 10 2 1 0 0
8 4 12.5 3 1 0 0
9 5 15 4 1 0 0
10 1 2.268 0 0 1 0
11 2 4.768 1 0 1 0
12 3 7.268 2 0 1 0
13 4 9.768 3 0 1 0
14 5 12.268 4 0 1 0
15 2 7.268 0 1 1 0
16 3 9.768 1 1 1 0
17 4 12.268 2 1 1 0
18 5 14.768 3 1 1 0
19 6 17.268 4 1 1 0
20 2 4.536 0 0 2 0
21 3 7.036 1 0 2 0
22 3 9.536 2 0 2 0
23 5 12.036 3 0 2 0
24 6 14.536 4 0 2 0
25 1 5.67 0 0 0 1
26 2 8.17 1 0 0 1
27 3 10.67 2 0 0 1
28 4 13.17 3 0 0 1
29 5 15.67 4 0 0 1
30 2 10.67 0 1 0 1
31 2 13.17 1 1 0 1
32 4 15.67 2 1 0 1
33 4 18.17 3 1 0 1
34 5 20.67 4 1 0 1
35 2 7.938 0 0 1 1
36 3 10.438 1 0 1 1
37 4 12.938 2 0 1 1
38 5 15.438 3 0 1 1
39 6 17.938 4 0 1 1
40 3 12.938 0 1 1 1
41 4 15.438 1 1 1 1
42 5 17.938 2 1 1 1
43 6 20.438 3 1 1 1
44 7 22.938 4 1 1 1
45 3 10.206 0 0 2 1
46 4 12.706 1 0 2 1
47 4 15.206 2 0 2 1
48 6 17.706 3 0 2 1
49 7 20.206 4 0 2 1
50 2 11.34 0 0 0 2
51 3 13.84 1 0 0 2
52 4 16.34 2 0 0 2
53 5 18.84 3 0 0 2
54 6 21.34 4 0 0 2
55 3 16.34 0 1 0 2
56 3 18.84 1 1 0 2
57 5 21.34 2 1 0 2
58 5 23.84 3 1 0 2
59 7 26.34 4 1 0 2
60 3 13.608 0 0 1 2
61 4 16.108 1 0 1 2
62 4 18.608 2 0 1 2
63 6 21.108 3 0 1 2
64 7 23.608 4 0 1 2
65 4 18.608 0 1 1 2
66 4 21.108 1 1 1 2
67 6 23.608 2 1 1 2
68 7 26.108 3 1 1 2
69 8 28.608 4 1 1 2
70 4 15.876 0 0 2 2
71 5 18.376 1 0 2 2
72 5 20.876 2 0 2 2
73 7 23.376 3 0 2 2
74 8 25.876 4 0 2 2
75 3 17.01 0 0 0 3
76 6 23.376 1 1 2 2
77 4 22.01 2 0 0 3
78 6 24.51 3 0 0 3
79 9 30.876 4 1 2 2
80 4 22.01 0 1 0 3
81 5 24.51 1 1 0 3
82 6 27.01 2 1 0 3
83 7 29.51 3 1 0 3
84 8 32.01 4 1 0 3
85 4 19.278 0 0 1 3
86 5 21.778 1 0 1 3
87 6 24.278 2 0 1 3
88 7 26.778 3 0 1 3
89 8 29.278 4 0 1 3
90 5 24.278 0 1 1 3
91 6 26.778 1 1 1 3
92 7 29.278 2 1 1 3
93 8 31.778 3 1 1 3
94 9 34.278 4 1 1 3
95 5 21.546 0 0 2 3
96 5 24.046 1 0 2 3
97 7 26.546 2 0 2 3
98 8 29.046 3 0 2 3
99 9 31.546 4 0 2 3
100 6 26.546 0 1 2 3


The answers are different in many instances. For example, to generate 34 cents, the minimum count solution uses 6 coins including a nickel, whereas the min-weight solution is 4 grams lighter and uses 4 pennies and 3 dimes.
MC:    34    6    20.67    4    1    0    1
MW:   34    7    16.804  4    0    3    0

Problem 2: Assuming that each desired value scenario 'v' is equally likely to occur, find an optimal distribution of coins to carry such that we can exactly satisfy each scenario.

We can model this by creating a vector of 'x' used for each desired scenario, as well as a single integer vector 'z' such that any 'x' value for a scenario is no more than its corresponding 'z'.

for(int i = 0; i < numCoinTypes;i++){
     cplex.addLe(0., cplex.diff(z[i], x[desiredValue][i]));
}

CPLEX log:
Tried aggregator 2 times.
MIP Presolve eliminated 45 rows and 41 columns.
MIP Presolve modified 4 coefficients.
Aggregator did 5 substitutions.
Reduced MIP has 450 rows, 358 columns, and 1067 nonzeros.
Reduced MIP has 40 binaries, 318 generals, 0 SOSs, and 0 indicators.
Probing fixed 0 vars, tightened 8 bounds.
Probing time =    0.00 sec.
Tried aggregator 1 time.
Presolve time =    0.00 sec.
Found feasible solution after 0.00 sec.  Objective = 395.3600
Probing time =    0.00 sec.
MIP emphasis: balance optimality and feasibility.
MIP search method: dynamic search.
Parallel mode: deterministic, using up to 4 threads.
Root relaxation solution time =    0.02 sec.

        Nodes                                         Cuts/
   Node  Left     Objective  IInf  Best Integer     Best Node    ItCnt     Gap

*    21     5      integral     0       36.5460       36.5460      969    0.00%

Thus, if we carry 10 coins (weighing ~36.5 grams) as distributed below:
penny:    4
nickel:    1
dime:      2
quarter:  3

we will be able to provide exact change for any scenario.

Problem 3: Determine the best 9, 8, 7, ..., coins to carry to minimize average expected absolute deviation from the desired values over all scenarios. Some of these turned out to be tough to solve to optimality due to the naive formulations employed. Nevertheless, optimal or near-optimal solutions were obtained in all instances.

Result: the optimal solution for n = 9, 8, 7, and 6 coins simply deletes a penny from (n+1) coin solution.
n =5, we delete a dime (0, 1, 1, 3)
n = 4, we delete a nickel (0, 0, 1, 3)
n = 3, 2, 1, are all-quarter solutions

Inventory Optimization
To more correctly model the residual change constraint (approximated via the objective function in Problem 3), suppose we are short by value 'delta': We can expend a dollar bill and end up with a net change of (100 - delta). In other words, we have to solve an inventory problem that for example, determines the optimal initial coin inventory vector 'z' such that the total weight of the expected final inventory after giving and/or receiving change over all scenarios  (and over multiple transactions or periods) is minimized. This model can be built by introducing a binary variable 'w' for each scenario to represent the case where a dollar bill is used or not used to supplement the value associated with 'z'. However, we do not know apriori, the distribution of the returned change, so some approximations are required. The analysis of this problem is a post for another day. 

Read part-2 here, and part-3 here.







Saturday, June 30, 2012

Alternative Optimal Solutions and Combinatorial Risk

Optimization in practice is usually not just about setting up a well-behaved model with an objective function and finding 'the optimal answer'. While that is an interesting exercise, the real 'value add' comes from the subsequent process of recognizing the business reality that a practical decision problem typically has many answers. Consequently, analyzing alternative optimal solutions in a way that is useful to the client is quite important. In other words, practical optimization is more often about analyzing feasible alternatives that initially appear to be equal to us, but in reality have vastly different qualities from our client's point of view (note: returning 'infeasible' is not really an option, and showing why our model returned 'infeasible' is only slightly more useful).  As we initiate a dialogue with our client to understand these differences and the context in which some alternatives are better than the others, we can see our lab model gradually transform into a useful business analytics tool.

Cutting across industries, I have not yet come across a single optimization problem deployed in production that does not have multiple objectives. Every seller loves to maximize profit, but not at the cost of losing out on volume or sales dollars in the process. Over time, the number and priority of such considerations change. For example, in the airline industry,  it is not uncommon for large-scale crew schedule planning problems to have hundreds of different goals and priorities. The richer the solution space, the more the number of goals it seems. In fact, optimizing just a single measure is risky because such a gain ("extreme point") almost always comes at the expense of other metrics that haven't been included within the analysis yet. Which leads us to:

Combinatorial Risk
This hidden problem of 'risk', even within a deterministic modeling context, is exacerbated in combinatorial (or global) optimization situations. Here, our model analyzes multiple inter-linked decisions that can produce solutions that radical differ from current practice and looks great numerically, but in reality, can potentially hurt our client's business if actually used. 'Locally optimal' does not always and automatically mean 'inefficient'. Like globalization, combinatorial or global optimization based holistic decision making can and does bring in more efficiency and profitability compared to that obtained by combining multiple locally optimal decisions when things go as per forecast. On the other hand, if the alternative (near-) optimal scenarios along with their corresponding risk of failure are not well mapped out and thought through, the resultant machine-generated combinatorial solution can cascade the risk of a bad decision through the system.


Part-3 here.

Thursday, June 28, 2012

Time-lapse view of a Soap Opera

One can speed-watch a bad 3-hour Shah Rukh Khan Hindi movie, set in some implausible post-modern world, in 15 minutes or less, and yet fully grasp the plot and story. Somewhere down the line, those wonderfully poetic song-and-artistic-dance filled world of authentically Indian movies and television programming disappeared, and along with it, my patience, leading to this compressed viewing method. On the other hand, it turned out to be quite interesting watching Jack Bauer save the U.S from annihilation in a '24' hour day, over a contiguous 24 hour period (on DVD) without any cuts. This provided a completely different experience and understanding compared to watching '24' live one-hour capsules across 24 weeks.

To experiment with a larger data set, I decided to watch a popular (5 days a week) Hindi TV serial that achieved high viewership ratings over its two-year lifetime (2008-10), and whose episodes were meticulously labeled and archived on the net by some fans. That data set had about 520 observations (episodes), which I speed-parsed over an exhausting six week period.

A TV Serial Storyline is a Random Process
Time-lapse viewing allows one to recall past events and key plot twists relatively more clearly, thus making it a little easier to figure out the 'rate' at which the storyline converges. If you think of a TV serial storyline as an origin-destination path that passes through a number of intermediate and interlinked stages (or states), you can identify the type of state transitions, including the many instances of 'back to square one' in a time-lapse mode. The soap storyline returns to a previously visited state (i.e. cycles or sub-tours) with a significant nonzero probability, much to the viewer's chagrin. You can also spot the random walks where the storyline drifts away from the original theme. When these events happen too often in soaps, it is a dead-giveaway that the plot is on life-support and the plug will be pulled shortly.

The story line of most TV serials appear to be remarkably Markovian in design. If you were to watch the series in real time (i.e. slowly), you would be able to skip many episodes and yet rejoin in a few weeks without really missing anything. In a soap opera, the future is indeed independent of the past, given the present.


Another O.R analogy
There are some intermediate episodes that generate peak ratings and a spike in viewer comments in the online archive. These often turn out to be events that appear to be 'pivotal' in the sense that they promise a lot of exciting ripple effects down the road. There is just the right amount of uncertainty in the future trajectory of the storyline but with some enticing structure around it. For example, if you were to look back at the Harry Potter book series, this would be around year 4 at Hogwarts. It's similar to what one experiences when solving a combinatorial MIP using branch-and-price (i.e. column generation) in practice. Starting with the root node, as you begin to apply restrictions and generate new columns, your LP objective function value usually drops*, and half-way through, your solver promises a really good solution. The possibilities are exciting. Alas, as your restrictions increase and you get close to the end state and your solution becomes increasingly integral, your objective function suffers a backlash from those nasty combinatorial constraints, and the objective function value goes through the roof. Invariably, much like a TV or book series, the ending is almost never quite as exciting as how one expects it would be. After investing all that time and effort. Perhaps, like life, good soaps are more about the journey and less about the destination.

*update (July 1, 2012)
The LP relaxation objective function value drops because:
a) the unrestricted subproblem in many large-scale planning problems is discrete, nonlinear, nonconvex; it is quite difficult to solve this problem to provable optimality within a reasonable amount of time

b) adding restrictions initially eliminates many bad solution areas, i.e. it simulates the role of inadvertent cutting planes, thereby improving solution quality.

Wednesday, May 23, 2012

Optimizing your Microwave Oven Performance using an Appalam

A far-out post on 'optimization' and 'parameter estimation'.

No two microwave ovens behave the same way even if they have the same power rating and capacity, and over time they show increased randomness in terms of energy usage and solution quality. Each one has its idiosyncrasies and as you move between apartments over the years, it is irritating to adjusted to a new oven. After you invest all the time in figuring out that the optimal settings to warm your cold coffee is 50 seconds in the previous home, using the same timing with the new one results in a lot of spilled coffee. Furthermore, objects get heated slightly quicker when placed in certain locations of the oven. What is a quick way to optimize your vessel placement and minimize your turn-around time? Enter, the Appalam.

The Appalam or the Pappad as it is called in Northern India is a circular, flattened, wafer thin dry mixture of lentils and spices. It is an almost fat-free, delicious snack if you microwave it although it can also be fried. In particular, the best brand for this test is the Lijjat Pappad. This one is hand made (or used to be). The Lijjat company rose from a tiny all-women cooperative start up in rural India (based on Mahatma Gandhi's principles of self-reliance, theirs is a remarkable and inspiring success story) and is still going strong, producing outstanding Pappads in a variety of flavors.
 

The idea is pretty simple: The Lijjat Pappad (LP) has a relatively large surface area that covers 75-100% of the circular glass tray in most residential microwaves. Microwave the LP for about 60 seconds and note which of its parts get heated up first (visually seen via a change in color) and how the cooking progresses. Within 60 seconds, you should be able to figure out the hottest and the coldest spots in the oven. In my current residence, the middle of the oven turned out to be the coldest, which was the exact opposite of the result in my prior residence. Of course, this is not intended to be a universal setting and only applies to a subset of foods being heated up.

Why the LP works relatively well for such a test:

1. The consistency of the LP mix is quite remarkable. It is neither thick to resist microwaving, nor too thin and very rarely exhibits any significant warping even after 120 seconds of microwaving (it will simply get carbonized before it warps).

2. The material naturally does not conduct heat well and convection doesn't help much either, so the localized heating effects show up visibly.

3. You can make a meal of your experimental subject once your test is complete. No test goes waste!

Wednesday, May 2, 2012

Optimal Shoelacing

There are gazillion alternative ways to tie your shoelaces of which only a few patterns are easy to remember and possess nice symmetry.


My daughter got her first laced shoes a few weeks ago. As she took it out of the box, the first noticeable thing was that the slack required to tie the laces appeared to be on the shorter side. This posed difficulties for my daughter while she was learning to make those neat knots. Rewiring the shoes resulted in quite a bit of slack, which resulted in big knots that not only helped her learn quickly, but also seemed irresistible to her kindergarten classmates who kept tugging at it. So the aim was to find a lacing pattern that would generate just the right amount of slack.


This is of course an instance of a Traveling Salesman Problem (TSP). Each lace hole represents a city which can be visited exactly once. For example, the slack is maximized by finding the shortest tour, while the slack can be minimized by finding the longest Hamiltonian circuit. Of course, my daughter would prefer finding a good balance between these two extreme solutions, while also ensuring that tightening and loosening the laces are relatively easy to perform. The former objective can be equivalently specified in terms of minimizing deviation from a desired tour length, while the latter requirement can perhaps be approximated by eliminating unfavorable connection patterns and reducing overall friction.


Any OR person will tell you this: just because the TSP is a notorious NP-Hard problem, it does not automatically mean that practical instances are terribly difficult to manage. On the contrary, OR methods excel in quickly finding amazingly (and provably) good answers to practical instances of underlying TSP substructures within decision problems across a variety of industries.

Saturday, April 28, 2012

House Hunting Inefficiently

If you are looking to be more time-efficient in hunting houses, see a prior post here

Halfway through the house hunt, the process began to get mechanical. There was not much diversity in house designs, especially the newer ones, although there was this really nice and shiny new house adjacent to a cemetery, which threw me off a bit. Internet forum opinions ranged from "Loved it. Quietest neighbors I ever had" to "Hey, its Halloween every night here!".  The older ones seemed to have more character and personality. However, they came with their own maintenance list. I then got quite interested in the principles of the amazing ancient Indian science of construction and planning 'Vaastu Shastra' that actually turned out to be a pretty useful practical guide (better than my realtor sometimes). And so the search continued ..

A common reason for many houses in the U.S (perhaps not as common in India) being put up for sale is that the children grow up and move away, and their parents want to downsize. Such houses are inevitably full of memory trails frozen within photo frames: childhood sketches, family reunions, high-school trophies, often ending with snaps of a daughter's wedding. Looking at those pictures changed the objective function. What was 'just' a house to tick off the list, was for many years a home where a family was raised from the cradle, parents aged gracefully, and kids grew up with security. That is no easy thing to pull off in today's world, and perhaps there would be very few things more satisfying that emulating what some parents in those families did. It was quite humbling. Each such 'ordinary' house had an extraordinary American tale to tell, some happy, some not so. And no matter how much money one pours into home improvements, that unique signature of how a family lived in that home does not really go away; after all, a family breathes life into a house. 'Must-have' product attributes no longer seemed that important. House-hunting stopped, and the search for a home began. This approach may be inefficient, but it certainly feels more rewarding and less tedious.

Saturday, April 21, 2012

Housing Industry Analytics: Short-Sale versus Foreclosure

There's been an explosion of 'short sales' of housing properties by banks in the US market. In the east coast markets, there's been a 60-100% spike in short sales, why? If the banks wait too long and then foreclose, the properties loses too much value and the final returns may be lower. On the other hand, a short sale is immediate, and long term risk is eliminated. As the article in the NECN link above says:

"In a "short sale," banks agree to let someone who owes more on their mortgage than their home is worth to sell it to a new buyer, with the bank typically writing off tens of thousands of dollars in the process. But what the bank gains is avoiding the cost, protracted process and uncertainty around taking the home or condo by foreclosure and then trying to resell it as a bank-owned property."

A first look indicates that this is a decision optimization problem under uncertainty that is somewhat similar to (but not the same as) that faced by fashion retailers who are trying to clear their end-of-season perishable inventory. Do they markdown apparel right now or should they wait for some more time? If they wait, then over time, the 'fashion statement' value deteriorates and the retailer may have to more aggressively markdown to attract customers and clear inventory. On the other hand, if they markdown right now, that may turn out to be a hasty and expensive decision, with a certain probability. So really there are two decisions to be made: when to markdown and by how much?

A bank may own a majority stake in several properties whose values are depreciating over time in an over-capacitated market (property owners are unlikely to have the cash to maintain or make improvements), which they need to get off their books without losing much. Using stochastic optimization methods available in the field of Operations Research ("the science of better"), they may be able to do a much better job of profitably managing their inventory (it won't be surprising if they are already doing this). Stochastic optimization methods are specifically designed to work with probabilities of scenarios, as opposed to a deterministic approach that assumes everything is perfectly known in advance, although the latter often turns out to be a reasonable and quick first approximation. An alternative that may be especially appealing to pessimistic banks looking to avoid worst-case meltdowns is 'Robust Optimization' that can operate without formal probability distributions and can among other things, help minimize a bank's maximum regret. Another advantage of OR methods is that they are typically not capital intensive, and the ROI on successful projects can be remarkably high. In short, OR can be very, very useful here.

Tuesday, April 17, 2012

Analytics and Cricket - VIII: DRS & Bayes Theorem

In the last post on cricket, we mentioned that the false positive (F+) issue with the Decision Review System (DRS) employed in international cricket could be a deal-killer (see red zone in picture below).

In this post, we work out an illustrative numerical example using a well-known conditional probability model based on reasonable data derived from interviews of ICC personnel to show that the current F+ rate disproportionally reduces the efficacy of the DRS, causing it to operate only marginally more effectively that the human-only (umpire) method, and thus may not be worth the cost of maintenance unless the F+ rate is reduced to a more acceptable level.

For brevity, let's focus on bowler reviews in this example. A bowler will ask for a machine review of an umpire's original decision of not out, hoping to turn that into an 'out'. Umpires in the elite panel are themselves around 90% effective in making the right decision (so on average, only 10% of the subsequent DRS referrals should change the outcome if they work perfectly), so it is really that 10% gap that is the problem.

Today's cricket DRS system is claimed to be around 95% accurate in giving a batsman out, if in fact, the batsman is really out. Suppose the DRS also yields F+ results for just 1% of the bowler reviews, i.e. it gives a batsman 'out' when he is really 'not out' just like the umpire originally said. If 10% of the batsmen subject to bowler reviews are actually out (as obtained in the previous paragraph), what is the probability that a batsman is actually out given that the DRS overturns the umpire's decision to say he is out?

Answer: Let OUT be the event that the batsman reviewed is actually out (its complementary event is NOTOUT), and RED the event that DRS gave him out. The desired probability P(OUT|RED) is obtained using the Bayes formula by:

P(OUT|RED) = P(OUTRED)/P(RED)
Expanding out the terms, we can write this as
= [P(RED|OUT) x P(OUT)] /
[P(RED|OUT) x P(OUT) + P(RED|NOTOUT) x P(NOTOUT)]

= [0.95 * 0.1] / [0.95*0.1 + 0.01 * 0.9]
= 0.095/0.104 = 91%

Observations
1. Even a 1% F+ rate brings down the true efficacy of DRS, and it is not 95% as the ICC claims. The second term is a combination of F+ rate and human accuracy. Thus
P(OUT|UMPIRE SAYS OUT) = 90%
If the bowler asks for DRS review:
P(OUT|DRS SAYS OUT) = 91%
Not much of an improvement

2. The better the umpires get at their job, the worse the existing DRS will statistically perform. For example, if the umpires improve their upon their accuracy by just one percentage point, i.e. to 91%, the conditional accuracy of DRS changes to:
= [0.95 * 0.09] / [0.95*0.09 + 0.01 * 0.91]
= 90%


Thus, the tables are turned now and the DRS makes things worse for batsmen here and we may be better off not using DRS at all even if it is provided free of cost!

This second result may seem puzzling. Why does this happen? If the umpires get better, the frequency of true NOTOUT is 1% higher, and with the F+ rate held constant at 1%, there will be an increase in the total number of false positives over a period of time, in addition to a small decrease in count of true positives, thereby reducing the accuracy rate of the DRS.

You can plug in a variety of numbers to see what the corresponding results are. You can also perform a similar analysis for batsman reviews.

Recommendations:
1. Significantly cut down on the F+ rate and not just focus purely on increasing true positive rate

2. Improve the quality of original human decisions. This will reduce the dependence on DRS, encourage improvements in the DRS to keep pace, and obviously improve player attitude toward umpires.

3. If a brilliant cricketing instinct filled person like Mahendra Singh Dhoni talks about 'adulteration of human and machine', do think twice about it, he's got a useful math model behind this statement!

Reference: Introduction to Probability Models by Sheldon M. Ross. This example is a variation of an example from this book. Hope I did not mangle it.

Monday, April 16, 2012

Anatomy of an Online Debate

The Huffington Post ran a curious online debate a few days ago: Is Yoga a Hindu Practice? Let me state right of the bat that I thought the debate was utterly idiotic given that this question was like asking people to vote if baseball was quintessentially American, if the great pyramids were Egyptian, or if the great wall was built by the Chinese. After all, how the heck does a person debate against a fact? By twisting it into a silly insinuation about ownership. Nevertheless, let's look at the debate results measured by market-share for, against, and neutral to the topic, tracked before and after reading the debate. 


If this picture is unclear:
For, against, neutral (baseline) = (65%, 26%, 9%).
For, against, neutral (after reading debate) = (66%, 28%, 6%)
So the net result is that the 3% of the undecided split 1% toward 'for' and 2% toward 'against' after reading the debate.

Warning: If you are a predictive analytics connoisseur or swear by rigorous statistical methods, the rest of the post will be cringe inducing, so read on at your own risk.

Given that we have absolutely no data to use other than this pie chart, I tried to 'quick fit' a plausible Multinomial Logit choice model to these results with the aim of personally understanding how useful this debate really was. Toward this, I defined a utility function u(t) = exp (a0(t) + a1(t)), where:
a0 = baseline contribution for (t = aye, nay, neutral)
a1 = debate contribution
market share (t) = u(t) / (u1 + u2 + u3)

Using the 'before' pie chart, I obtained the following values:


where CONST = a0, and a1 = 0 at this point. Again, note that these are just plausible values and not statistically calibrated likelihood maximizing coefficients based on historical individual observations. Next, merrily using the  'after' pie chart to update the utility function taking the debate into account, I obtained the following values:



where DEBATE = a1 that is used to update the utility function. Note that the results don't really change dramatically.

Based on this plausible MNL model, we observe a positive value for a1 for both the 'for' and 'against' since their market-shares increase after the debate, and a negative value for 'neutral' since the debate forces a good chunk of the few fence sitters to switch. To measure the usefulness of the debate to each group, I looked at the ratio of "what was additionally useful" versus "prior understanding", i.e. the ratio utility_before/utility_after given in the last column. The results indicate that the debate itself was pretty close to useless to the overwhelming majority of the voters, i.e. for the 'aye' people (like me) and on the whole, the debate reinforced what they already knew, yielding a tiny usefulness change of about 0.4%. On the other hand, the debate was relatively more useful to the 'nay' people and it incrementally 'hardened their position' by about 6%. More than a third of the miniscule fence sitters actually took a stance and this debate appealed most to them.

Given that the voter response ('elasticity') for a particular choice-group to an event (debate) in an MNL model also depends on the incremental gain possible from their existing market-share, the degree of movement in market-shares for the three groups are not surprising, although the specific direction of the resultant net shift in market shares does indicate that the debate may have had something to do with it. Of course, one can game such online "changing minds" debates by entirely ignoring the debate and starting with a non-favorite position as the baseline and then simply selecting your most favorite position in the end.

On a side note, almost all the 'Yoga' practiced in the U.S and the west is really Yogasana, whose primary function is to help prep your mind and body for actual Yoga, which in turn has nothing to do with whether you can twist yourself into an exclusive USPO patent-protected double-pretzel or not, and everything to do with open-source inner-sciences that aim to rid a mind of ego, exclusivity, and dogma and reach higher levels of consciousness.

Monday, March 26, 2012

Gender-Shaping III: Is Amartya Sen's Missing Women Count Exaggerated?


This is the third post in this series on gender-shaping. The previous installment can be found here. Thanks to a twitter link, I came across a 2010 journal paper: "Missing Women: Age and Disease," Siwan Anderson (University of British Columbia) and Debraj Ray (New York University) published in Review of Economic Studies Vol.77.

This paper has among other things, investigated Amartya Sen's '100 Million Missing Women of India' claim that is attributed to systemic discrimination. Anderson and Ray have estimated the number of 'missing women' in India, China and Sub-Sahara Africa by age and cause-of-death (not done before) while also moving away from the simplistic aggregate sex ratios that were used as baselines in prior works. The authors make the following useful observations: Defining missing women by differences in aggregate sex ratios can be misleading, or uninformative (or both). It is misleading because different countries have different fertility and death rates, and (in particular) different age distributions. They will have different disease compositions.
They may also have different sex ratios at birth for genetic or environmental reasons that have nothing to do with missing females
.

The procedure is also uninformative: we cannot tell at what ages the missing women are clustered, or what diseases are responsible. Thus, we cannot begin to ask about the various
channels: discrimination, biology, social norms, and so on. Answering these questions is of profound importance. By unpacking missing women by age and disease, our paper takes a limited and preliminary step in this direction.
"

From an OR perspective, we extensively rely on similar customer segmentation models (in revenue management for e.g), and this additional age- and causal-factor based segmentation appears to be quite important and yields two main results as well as a comparative result that may be interesting to an U.S audience:

1. A large fraction of the missing women in India are not infants (less than 20%) but adults, and is attributable to other factors like disease and injury, apart from any systemic discrimination. Consequently, any claim of exclusively female infanticide driven 'missing women' in India is rejected. On the other hand, this paper finds that 44% of China's missing women are in the prenatal age-group. Here is a snapshot of sex-ratio by age, taken from the Anderson & Ray paper:




2.The authors make an interesting comparative comparison with the U.S: "we observe some similarities between age-specific percentages of missing women in the historical United States (ca. 1900) and India or sub-Saharan Africa today".

3. The Sen count (100 million missing women) appears to have been calculated with respect to a specific counterfactual: The overall sex ratio for N. America, U.S and Japan. An alternative calculation by Coale (1991) comes up with a more conservative estimate of 60 million. Anderson and Ray perform similar calculations but at the segment level (i.e. by age-disease) and generate missing number estimates using more carefully chosen counterfactuals as the baseline and find approximately 20 million missing women in India (aggregated across all age groups), while the corresponding figure for China is 58 million. Furthermore, 'injury' is not an insignificant culprit in India across all age groups, a potentially worrying trend that its government must look into. (The paper alludes to the old bogey of 'dowry deaths' as a probable cause which may not turn out to be the case. A similar detailed analysis is required).

The findings of this paper also weakens a statement in a previous post on this topic that a skewed overall male:female ratio in a region is a 'scary indicator' of female infanticide being practiced there. My statement ignored the age distribution as well as the 'cause of death' dimension. Bad O.R, but I have Amartya Sen for company.

Thursday, March 22, 2012

The Optimal Playlist

One of the problems with neighborhoods in parts of Connecticut is the lack of sidewalks coupled with crazy drivers (probably from a neighboring state to its left). To avoid getting run-over, I decided it is safer to do my walking on the treadmill. I'm now getting all the exercise a creaky researcher needs, but I'm not getting anywhere. To overcome this monotony, I hooked up my old iPod-classic for company, but it's time-consuming to generate my preferred playlist : start off with some up-temp music for motivation, then switch to cruise mode, and tone down after my 30/60 minutes of walking.

In India, we have the concept of 'Rasa', a Sanskrit untranslatable that very roughly speaking, includes notions of experiencing certain emotion(s), themes, ambiance, genre, etc. So the sequence of Rasas  matters a great deal. Furthermore, I like to listen to complete songs and hate to end a virtuoso Carnatic performance half way when the exercise session-clock runs out. Furthermore, there are so many languages in India and many have their own pop-culture, folk, and classical genres in instrumental as well as vocal modes, and I prefer a diverse sampling of these to feel more at home.

Putting all this together to achieve an optimized playlist requires a constraint-programming approach. If I also want to optimize a certain objective (e.g., stay close to 12 songs), this turns into an exercise of solving an associated discrete decision optimization problem that can be stated as follows:

Find (preferably) 12 complete non-repetitive songs in a preferred sequence that lasts (almost) exactly 60 minutes, and includes at least n(i) songs having user-specified attribute (i), i = 1, .., n.

If we restate the attribute requirement as a soft-constraint by creating a score-table for including any attribute (e.g. 10 points for including a song with attribute i once,  15 points for two songs with attribute (i), 17 if three or more times) as opposed to the 'must satisfy' version stated earlier, then the playlist optimization problem can be posed as a attribute score-maximizing multiple-choice knapsack problem with a cardinality constraint, followed by a sequencing step. Even with a huge home music database, practical instances of the latter formulation may be relatively easy to solve via combinatorial methods (iPhone app?) and may not require expensive MIP solvers. Then, as a second step, we can sort the included songs into another preference-score maximizing sequence to generate the final playlist, unless of course the sequencing requirements are not that simple (in which case, a more sophisticated optimization approach may be required).


Such an optimized playlist is also useful if you want to build an auto-pilot DJ for your next house party. If your approach can solve this problem on-demand, you would also be able to dynamically re-optimize the playlist after manual intervention.

It seems apt to terminate this post with a Carnatic-Western classical fusion piece.



Updated on March 30: The objective function above is deterministic so there is a good chance that the you will get the same set of songs to listen to each day, which is not very useful. To introduce diversity and exploit the fact that in practice you tend to get several alternative optimal solutions to such problems, add a small amount of clock-dependent noise to the attribute-score and sequence-preference score. This will likely do the trick.

Thursday, March 1, 2012

House Hunting Efficiently

One of the consequences of the kind of convergence shown in this graph was that it created the need to buy a house. It's become a ritual to spend weekdays creating a list of houses that are feasible with respect to hard constraints (big kitchen, level lot, ..), and then converting that into a prioritized list based on how they score on soft constraints (pre-wired for Bose speakers, for example) that in turn motivates a preferred way of visiting these houses during weekends. I noticed that I rarely see each house in isolation and my view of a house tends to be colored by what I saw earlier. However, of late, the time to view these houses has become a scarce resource, so I created an O.R driven prioritized list that maximizes and optimally allocates viewing time, keeping the total duration equal to the limited time available. I used my automobile GPS unit as the solver.

This GPS unit "solves" the Traveling Salesman Problem (TSP) to figure out the optimal (?) order of visitation that minimizes total drive time, which automatically maximizes aggregate viewing time. (In particular, if houses are located on either side of a busy highway or Main Street, a good heuristic would be ensure that the optimal path intersects such a link infrequently.) I can then allocate the optimal expected viewing time to the houses based on personal preferences. The total viewing time also informs me if I have spread myself too thin, in which case, I can start deleting houses with the lowest scores from the list and re-optimize until the solution looks reasonable. 

If viewing order is important from a 'relative comparison' perspective, the resultant constrained TSP problem becomes a bit more harder to solve using a GPS unit. A simple heuristic rule could be to fix the second node ("first house to first") and/or the second-last node ("house to visit last") of the tour and let the others be visited based on time-optimality. If your realtor drives you around, her/his office is the start and end node of the Hamiltonian circuit.

One issue I encountered while using the GPS unit to merely drive-by a house as part of a local neighborhood search (pun unintended) is that I had to get close enough to the house and maybe pause a bit to inform the GPS that this house has been reached. Otherwise, the GPS unit would continually re-route me back to the house, resulting in considerable confusion.