Cfrgtkky

Say I have a matrix

x1 x2

x3 x4

With x1, x2, x3 and x4 randomly and uniformly drawn from the interval [0,1]

After I do gauss-jordan elimination, what is the expected value of the absolute value of the second pivot (the one on the x4 position)?

So far I seem to be getting a surprinsing answer: infinity

For the second pivot can be calculated as $x4-(x3/x1)*x2$. Therefore, the sought expected value is $E(x4-(x3/x1)*x2) = E(x4)-E(x2)*E(x3)*E(1/x1) = 0.5-0.5*0.5*E(1/x1)$.

But $E(1/x1)=infty$ (see Expectation of 1/x, x uniform from 0 to 1).

Is there an error somewhere in my thinking? Where can I read a bit more about this?

Your computation appears to be correct to me. It's also supported by the following little bit of matlab code:

function gjtest(n)

mats = rand(4, n);



cells = zeros(n, 1);

for i = 1:n

    s = reshape(mats(:,i), 2, 2);

    s(2,:) = s(2,1)/s(1,1) * s(2,:);

    cells(i) = s(2,2);

end

mean(cells)

clf;

plot(cells); 

figure(gcf);

which produces, as mean values:

n       mean

1000     2.17...

10000    2.14...

100000   2.94...

1000000  3.100...

10000000 4.441...

which is not growing very fast, I admit, but certainly suggests that the mean is not $1$.

So you have a random variable with a very high mean, but the odd characteristic that its "population mean" tends to be a great deal smaller. I think this indicates lots of "skew" in the data (I apologize if I'm misusing statistics terms...it's been a long time). The plots from running those code fragments certainly suggest a very biased distribution. If you replace "plot" with "histogram" you can see what I mean.

In fact, let me add a little more. The pivot $X$ is (up to small perturbation), roughly $1/x_1$. How is this distributed? Well,

begin{align}
P(a-u < X < a+u)
&= P(a-u < frac{1}{x_1} < a + u)\
&= P(frac{1}{a+u} < x_1 < frac{1}{a - u}\
&= frac{1}{a-u} - frac{1}{a + u}\
&= frac{(a+u) - (a-u)}{(a-u)(a+u)}\
&= frac{2u}{a^2-u^2}\
&approx frac{2u}{a^2}\
end{align}
for small values of $u$. The pdf, gotten by taking a limit of difference quotients, therefore gives
$$
p(x) = frac{1}{a^2}
$$
This explains why there are so few large values of $X$ -- the probability falls off quadratically with the target values. But it also explains why the expected value is infinite: if you look at $x p(x)$, whose integral is the mean of $X$, you find you're integrating $frac{1}{x}$, which gives you a $log$, and the integral diverges. Then again, it diverges about as slowly as something possibly can [I'm speaking informally here!], so it's no surprise that the values I was getting were in the $2$ to $5$ range. :)