Topology?

Topology is not my best math.  Sure, I solved the bridges of Koningsberg problem in high school as you do, and I can feign a laugh at the doughnut/coffee cup trick that isn't decades old and seriously come up with a new joke already, and I realized what the shape of the Space War universe was while playing the game on an Atari 2600 with no sun and no bouncing (also: old video games always had kick ass covers).

Why are we talking about topology?  In E17, I like to have a number of virtual desktops, and I like them to be all connected to each other so I can flip between them by pushing on the edge of the screen.  FVWM did this years ago, and I always liked the idea.  E17 also supports wrapping of virtual desktops, so you end up with something like this:

My laptop uses a three by three grid, but it's the same thing.

Where there is a two-by-two grid of desktops, each connected to the others along wrapped edges.  This is obviously identical to the SpaceWar universe, as you can imagine this as just a really crappy tv with only four pixels.  Therefore, my standard desktop arrangement is homeomorphic with a torus.  Easy peasy.

Now, at work, I have a dual monitor setup, and I assumed this would just work the same, but with fat desktops.  However, I discovered upon getting everything configured that it's actually this far more complicated thing:

Green arrows are bidirectional, as they were above.  Orange are unidirectional.

The two monitors have their own set of virtual desktops, but, since going across the monitor boundary must put you on the other monitor, they're not connected together.  If you're on the left monitor, there is no right edge you can cross to get to left side of the adjacent desktop (imagine the center red screen to pink).  You can get there by going to the unopposed left edge of this desktop to arrive on the right edge of that desktop, but you then can't go back, as you jump over to the other monitor.

I don't have proof (again, see "worst math"), but I believe this defines a very oddly connected hollow torus.  I don't know how you construct a unidirectional connection in topology, but assuming that's valid, you basically have some sorted of ratcheted torus for each monitor: you can travel poloidally unobstructed (corresponding to vertical shifts), but can only travel in one direction toroidally.  Switching between monitors is then functionally equivalent to passing through the ratchet onto the other side of the torus.  There's some ambiguity, since you can then shift each monitor separately, so where you come out on the other side of the torus isn't a priori obvious, but I think the bulk of the thing is sorted.

Anyway, I thought it was cool, and now my already-confusing-to-other-people desktop is going to get more confusing.

I am irrationally concerned with good statistics

k, statistics again.  The problem is that I saw this article today, which basically complains that "no one really means to use standard deviation, as people intrinsically want to use the mean absolute deviation" which is, of course, completely dumb.

First, no one would ever do mean absolute deviation in their head.  Here are some numbers: {-1 2 3 -5 1 400}.  If you had to guess another number that would belong to this set, you're going to guess like "dunno, zero maybe?"  You know that 400 is probably wrong, so you cut it out.  People don't do real means when they filter data.  It's some combination of a mode and median.  Choose a number that doesn't seem crazy.

Second, this mean absolute deviation tells you about where the 50% point falls.  Why that point?  The standard deviation is more inclusive, as it tells you that most (Q(1) = 68.change%) samples are closer to the central value.

Third, all that obvious stuff about moments analysis.

Anyway, time for plots.  These are the same idea as the ones from the previous post, just remade with more samples and different stats.  The horizontal lines are the true uncontaminated distribution sigma and the true fully contaminated sigma (sigma_uniform = sqrt((b - a)^2 / 12), because math).  First thing to note:  Actual sigma cleanly switches from the two extremes, as it really should.  Gaussian fits are best, but IQD and MAD are comparable up to the 50% contamination point.  MeanAD doesn't seem particularly good.  The full contamination end is biased, as I'm using a parametric model (that it's a Gaussian distribution).

Biased samples.  This nicely shows that IQD fails before MAD, and that Gaussian fits are reasonable up to 60% contamination.  MeanAD is again off kind of doing its own thing.  Median >>> mean for outlier rejection.

Fish sandwich is about this early, and twitter works right on tv because I see how the lazy weekend is bizarre conniving prostitutes prostitutes

-Or-

So, let's talk about Markov Chains, I guess.

This news story popped up in the RSS today, which reminded me of the long ago time in which Forum 2000 was a thing (note: Forum 2000 never worked like that.  It was all basically a mechanical turk).

I then remembered that you can download all your tweets from twitter, so blammo, I have the text corpus necessary to hack up something to do this.

So, those Markov Chains.  Here's the simple way to think of it:  starting from initial point A, there's a probability p_B that you'll move to point B, p_C that you'll move to C, etc.  However, there are also probabilities for all those points too, so you can chain things together.  Maybe go look at the wikipedia figure.  This isn't turning out to be that simple of an explanation.

In any case, if you let each "point" be a word in a tweet, then if you have a large sample of tweets from someone, you could imagine that you could construct fake tweets if you knew the probability that a given word follows another.  That's what I did in perl, and that's where the title came from.  I did a simple chain, where w_{i+1} is drawn directly from the probability distribution of words that follow word w_i.  There's a reset condition where if the set of words that follow w_i is empty, I restart the chain from the set of "words that start a tweet".  I've also forced the length of the generated string to be 140 characters, the standard tweet length (and forced everything lowercase, and removed all punctuation, etc, etc, etc).

Here are some more examples:

• i was far the whole dinner but that @wholefoods sushi was aborted because its weird being even fucking clue what hes a cheater :-p theres no
• hey verizon network that too hungry or more not sure that flasks were afraid the sun orbit simulator its now its someplace else too soggy as
• happy birthday @jkru isnt helping anybody know and hes never a 12 16 garbage bags tied directly to work - british accent for being wasted its
• rt @vgc_scott man that quote and cranky tonights dinner oh man its still preparing the oil is actively looking up a workaround to be going to
• first time you dont want a walk an impenetrable island which is like it looks close eyes will always just woke up dry before i not having intermotrons
• @jkru its all the music #fuckedupdream peter potomus cartoons are you your app its official word for liberty tsa approved for my implementation
• dear this is what is terrible monsters outside who like fog and colby isnt even left there to the 100th day yay deposed monarchy ok so hopefully

Those are just nearly not gibberish enough for you to think that they're not computer generated.  The crazy thing about MCs is that they do a wonderful job of constructing sentences that are nearly grammatically correct.  This could probably be improved a lot, such as trying to pick words that have some influence from the second previous word or having it clearly mark where it's restarted the chain with a period or something.  Still, for like ten minutes of work, this isn't too bad.

Finally, an interesting sampling from printing out the chain transition probabilities:
fuck you 0.294392523364486 63 214
fuck it 0.0794392523364486 17 214
fuck fuck 0.0700934579439252 15 214
fuck yeah 0.0700934579439252 15 214
fuck up 0.0327102803738318 7 214

[...]

fuck #sandwichtweets 0.00467289719626168 1 214

fuck shitshitshitshitshitshitshitshitshitshitshitfuckshitshitfuck 0.00467289719626168 1 214

[...]

fucking ghost 0.00454545454545455 1 220

[...]

fuckity fuck 1 3 3

Histogram binsize

Given a list of values you want to take a histogram of, what is the best binsize to use?  For data that's basically flat, you can just set the end points at the minimum and maximum values, and choose how many points you want per bin.  Assuming that the noise per bin is basically Poissonian, you can define the S/N you want for the bins, and then note that you then want on average k points per bin:

S/N = sqrt(k)

That's the easy case.  How do you do the same thing if the data is drawn from a Gaussian?  It's basically the same kind of an issue, except you need to choose where you want to optimize the S/N.  For the case from work, we largely ignore everything outside of 2\sigma, so I've used that as the critical point.  Now, for N total data points, how large of a box do you need centered on the critical point to achieve that S/N value?  Conveniently, this is just an exercise in error functions:

k/N = normcdf(critical + binsize/2.0) - normcdf(critical - binsize/2.0)

k/N = 0.5 * (erf(sqrt(2.0) * (critical + binsize/2.0)) - erf(sqrt(2.0) * (critical - binsize/2.0)))

So, you choose the S/N value for the critical point, determine the k/N value, and then find the binsize that achieves that (such as via the following figure).

This is for critical = 2\sigma.

Linear algebra problem.

It's not immediately obvious, but any sort of solver, fitting, or minimization is all basically just linear algebra.  This week's work-related math problem was in the polynomial fitting code.  It only arises when the x-data is closely (like 1e-3) spaced.  How is this linear algebra?

Here's how you fit a polynomial.  First, define the Vandermonde matrix for a vector x as

VM[x] = [x.**n x.**(n-1) ... x.**2 x.**1 x.**0]

using matlab/octave style .-operators for element by element operations.  Basically, this has expanded the vector x into an array that has the various power terms in the columns.  You can then define a system of polynomials as VM[x] p, where the vector p contains the polynomial coefficients.  In the fitting form, this vector is what you're trying to find.  To do this, you do basically the same thing as regular fitting:

A p = b

(VM[x]^T VM[x]) p = (VM[x]^T y)

using some prior knowledge of how least squares fitting usually goes.  This is basically a matrix inversion (of A), but that's hard to do, so usually this kind of thing is solved with a LU-decomposition or the Gauss-Jordan method.  However, as mentioned above, if the x-data is closely spaced, this is failing, as our GJ implementation has some logic to prevent the construction row operations that have large scaling factors.

This leads to the idea of pre-scaling the data in x to prevent this from being a problem.  Define a new vector

x' = f x

where f is some scalar.  This means that

VM[x'] = VM[x] .* VM[f]

so

A' = (VM[x]^T VM[x]) .* (VM[f]^T VM[f])/length(x)

b' = (VM[x] .* VM[f])^T y

which has the solution

p = p' .* sqrt(diag(VM[f]^T VM[f])/length(x))

Great.  What's f?  This is where things start to be complicated.  Why does the GJ solver fail?  It's complicated, but basically boils down to the matrix A being ill-conditioned.  The degree of ill-conditioning for a matrix can be calculated as the condition number (with the most ill-conditioned matrices being non-invertable, where c_N = infinity).  The fastest way to do this is to calculate

c_N' = max(diag(A))/min(diag(A))

If this is large relative to some threshold (like 1e7), the matrix is ill-conditioned and should have some scaling.  Unfortunately, this is just a decent lower limit, so to do this correctly, you'd need to do something like

c_N = abs(max(eig(A))/min(eig(A))

or, even better:

[U S V] = svd(A)

c_N = max(diag(S))/min(diag(S))

In any case, I have the suspicion that the "optimum" f value is something like

f = sqrt(max( max(diag(S)), 1 / min(diag(S))

The idea being that you want to reduce the dynamic range of A in such a way that you no longer break the threshold.  One complication of this is that checking on the real matrix that's breaking things suggests that the "scaled" form that works is actually even more ill-conditioned than the unscaled, but since the test in the current code is uni-directional, this doesn't throw an error.

That was kind of hard to go through, so here's a Jigglypuff mocking Clefairies.

Possible solution to a major problem.

Allowing me to do this (left: original; right: modified reconstruction).

Given image i, such that i = a + n + c, where a is the astronomical signal, n is the random noise, and c is the detector specific corruption, define a new image i' = (a - a') + n + c, where we've removed an estimate of the astronomical signal a'.  This has a Fourier transform F(i') = (A - A') + N + C, using the fact that Fourier transforms are linear.  Assuming a' is a sufficiently suitable proxy for a, and that n is small (specifically that ||N|| << ||F(i')||), this leaves F(i') ~ C.  Constructing a Fourier mask that removes C results in t ~ (a - a') + n + c.  The difference of i' and t provides a clean estimate of c, which can then be removed from i: o = i - (i' - t).

I need to test that photometric properties are retained by this transformation, and I'd like to speed it up some and sort out some of the edge case issues (visible around the blank bands above), but I think this solves the bulk of the issues with this problem.

The Supreme Court

Specifically, killing them.

Hang on there, Secret Service.  We're doing it statistically, so it's just an math problem.  No one is actually getting hurt, and I mostly did it to see how soon it's likely that there will be a vacancy due to death.  Someone could resign at any time, and that would make a vacancy, but that's not something that can really be predicted.  Death is a bit tricky, too, but due to the fact that despite the best attempts of some people, we still have a government that collects piles of data, so we can use the cohort life expectancy data from the Social Security Administration.

From that data, we're going to use the "High-cost" columns, and use the "at age 65" data.  The logic behind this is based on two facts.  First, the majority of the difference between the "at birth" and "at age 65" expectancies is due to infant mortality and things like that.  If you survive to 65, you didn't die as a baby.  The supreme court is universally comprised of people who did not die as a baby.  The choice of the "High-cost" is that being a supreme court justice is a pretty good job, so they're likely to have good health care.  The final assumption I have to make is that people born before 1940 have life expectancies similar to those born in 1940.  The table only goes back so far, and I don't really want to extrapolate it.

The next assumption is that the probability of surviving a year is based on the logistic function, largely due to a lack of better options.  The form I chose was P(age) = 1 / (1 + exp(-(expectancy - age))).  This sets the 50% probability at the expectancy age.  I could scale it, but I don't have any good idea of how to do that, so I'm just letting that be 1.0.

Using this table of birth dates and ages, we can make the following plot:

This shows the probability that that justice is alive.  I'm sorry that this claims Ruth Bader Ginsburg is 52% dead.

Ok, using these probabilities, we can plot up the number of justices statistically alive by simply summing these curves:

This suggests that (barring the Senate being jerks), Obama is likely to appoint two justices.  His successor falls into a generational gap between the Breyer and Thomas, and isn't likely to appoint anyone, even in an eight year term.  It's fairly linear over the next five terms, with about one each term.  This assumes that all newly appointed justices are younger than the current youngest, and that we're only concerned with the current set.

So again, barring the Senate terrible, we should have a significantly more liberal court within the next four years.

Batman and Robin are happy about that fact, too.