Review of natural field transfer function algebra

Learning Objectives:¶

• the basic TF Equation set-up and how OLS solutions solve it
• How to cast the solution in terms of cross-and-auto-powers of array channels
• The concept of MT signal and noise
  • the difference between coherent and incoherent noise
• The role of incoherent noise in results bias
• The role of the remote reference station in reducing bias

What is not covered:¶

• We start with the concept of an MT Station and time series of calibrated Fourier coefficients
  • No discussion in this module is provided for how these are obtained

In a noise free world, at a fixed frequency, under the plane wave assumption, the measured fields at any fixed location should obey the following equation:

The fundamental TF relationship is defined by equation \eqref{eq:1}. The elements of the 3x2 tensor $Z_{ij}$ and $T_i$ are (theoretically) invariants for a fixed measurement location for simultaneous measurements of $E.

*Note that we can regress on $H_x, H_y$ as if they were “inputs” and the $E_x, E_y, H_z$ were “output” even though physically the $H_x, H_y$ may be some combination of input and output.

This is three equations in six unknowns, but each of the three equations can be solved separately for pairs of tensor elements (rows of Equation 1).

Extracting the first equation from the system we write:

In the “E_x” equation. the values for E, H are themselves vector observations

$

$

This is an equation with two unknowns, $Z_{xy}, Z_{xy}$.

Here the data are organized into row vectors. $E_x$ is 1 x $N$, $Z$ is 1 x 2, and $H$ is 2 x $N$. ¹

The standard (Ordinary Least Squares -- OLS) way to solve this is by multiplying $H$ on the right by $H^{\dagger}$ , the $N$ x 2 conjugate transpose, yielding a square² matrix, and then multiplying on the right by the inverse of this newly formed square matrix.

The next few equations parameterize the OLS solution in terms of cross powers, but they boil down to the following, 2d solves 2c for OLS:

$

$

$

$

Footnotes¶

¹ If we used column vectors the equation would be modified so that Z is right multiplying H.

² However, in practice, all that is required to solve the equation is to muliply by a matrix C such HC is invertible.

The OLS Solution¶

The matrix multiplications in \eqref{eq:EXTFHDaggerSolutionStep} involving $H^{\dagger}$ can be written in terms of cross powers

Let $\left<A B^* \right>$ denote the cross power (inner product) between complex vectors $A$, and $B$

Then we can express the above as:

The equation is now in a form where we can call numpy or scipy’s solve functions, but note that since the solution is only $2x2$ we can also solve explicitly.

The 2x2 on the right (we can call it $HH^{\dagger}$)is basically the covariance matrix of the magnetic field data.

A has the following inverse:

where the determinant is:

Substituting the form of 4a into 3b we have:

Now multiply both sides by 4b (using the expanded form on the left)

On the right hand side we have isolated our solution as $HH^{\dagger}$ cancels with with it’s inverse (assuming it exists!)

Expanding the LHS by multiplying out the cross powers we arrive at the solution:

We can now express separate equations for $Z_{xx}$ and $Z_{xy}$

And, finally, substituting 4c for the determinant:

Now we substitute 4c into the denominators to get $Z_{xy}$:

And for $Z_{yx}$:

N.B. We recover the solutions for $Z_{yx}$ and $Z_{yy}$ simply by substition of $E_y$ for $E_x$ in the above equations, and similarly $T_{x}$ and $T_{y}$ by substition of $H_z$ for $E_x$. Proof is left as an exercise for the reader.

A more general approach¶

In Equation 3a, the standard method of using $H^{\dagger}$ was applied, but what happens if we use some arbitrary matrix $C$,

Go back to Equation 3a, but instead of using $H^{\dagger}$, use $C$, made of column vectors A, B

The solution proceeds identically, (c.f. 9a with 3a)

But now the analog to 3b is :

Now the 2x2 on the right is the covaraince matrix between $H$ and $C$

With corresponding inverse:

and determinant:

Substituting the form of 10a into 9b we have (cf 5a):

Now multiply both sides of 11a by 10b (using the expanded form on the left, and the condensed on the right) to get Equations 11b (c.f. 5b)

Again isolate our solution as $HC$ cancels with its inverse (assuming it exists!!!)

Again Expanding the LHS by multiplying out the cross powers

We can now express separate equations for $Z_{xx}$ and $Z_{xy}$

Finally, expanding the denominator, we have separate equations for $Z_{xx}$ and $Z_{xy}$ in terms of cross-powers

Summary of the Single Station Formualtion¶

We have shown that the elements of the impedance tensor can be derived from the cross-powers of the various channels.

Note that this has all been one under the assumption that the measurements are being made at a single location (that is what the TF equation is relating -- Ex, Ey, Hz, to local Hx, Hy.

Mathematically, any channels $A$, $B$ can be used to invert the matrix, provided that $HC$ is invertable, which is the same as saying it’s determinant is non-zero. For example, we cannot allow A equal B.

• TODO: There is a stronger condition which also must be obeyed, which is that the columns of $C$ must be coherent with H, E.

Sims et al. evaluate each of the 6 pairings of Ex, Ey, Hx, Hy (from a single station).

• They point out that pairings (Ex, Hx) and (Ey, Hy) are unstable -- this is because tend to be co-linear, so the determinant of $HC$ can be very close to 0.

The remaining combinations are all “legal”, with two of them tending to overestimate and two tending to underestimate the impedance tensor elements, ... this can be seen by auto-power and cross power arguments (TODO add more detail here -- AST).

A corralary to this that seems worth mentioning is that any rotated version of the channels $A$, $B$ would still result in an invertable HC. i.e. $C$ is Nx2, so we can muliply it by any 2x2 rotation matrix and we still get a linear-independent Nx2. Say that we represent the rotation matrix by $R$, then $(HCR)^{-1}$ = $((HC)(R))^{-1}$ = (R)^{-1}(HC)^{-1}$, i.e. we will recover the same solution no matter if the reference channels are rotated.

Noise and Bias¶

It will help to consider a 1d earth and perfectly orthogonal sensors... in that case $Z_{xx}$=0.

• Proof left as an exercise for the reader -- Hint: use Assumption 0 below.

This is basically equivalent to Assumption 0: $H_x$, $H_y$ are independent and incoherent.

• For a homogeneous earth this implies (by the Maxwell-Faraday Equation) Corollary 0 that parallel $E$ and $H$ are incoherent.
  • Exercise: Does this also hold for a 1-D earth? How about general 2-D? What about 2D along-strike-rotated?

So the first set of equations reduces to $Z_{xy}$ only, .. repeating Equation 13b here

If we used the standard local channels, $A$ = $H_y$, $B$ = $H_y$, Then $Z_{xy}$ becomes (repeating equation 7b):

But by Corollary 0 the first term in the numerator is zero, (so is the second term in the numerator by Assumption 0 -- Hx and Hy are independent), so:

by Assumption 0 first two terms in the denominator are also zero, so:

Cancelling minus signs and the autopower in Hx

Equation tmp2 will be used to illustrates why some channel pairings (Hx, Hy) resulted in an under estimation of $Z_{xy}$ ala Sims et al. 1971.

Noise¶

What if the data are not noise-free? What can we say about the derived solutions?

That would depend on the details of the noise, so start with a simple noise model based on a few assumptions:

Let any channel of data, $a_{obs}(t)$ = $a_{s}(t) + a_{n}(t)$. ³

The subsripts $s, n$ are for the “signal part” and the “noise part” of the observed data respectively.

Let’s evaluate the role of noise under some Assumptions:

• Assumption 1: Noise is incoherent from all signal
  • $<a_{sig}(t),b_{n}(t)> = 0$ for any channels $a$, $b$. We haven’t really formlalized what is signal and what is noise yet, so this is a pretty weak constraint. For now say that signal plane wave, that obeys equation 1. So all noise is assumed incoherent with the signal source.
• Assumption 2: Noise is incoherent between channels
  • $<a_{n}(t),b_{n}(t)> = 0$ for any channels $a$, $b$. So the non-plane-wave part of the observed data that does not obey equation 1 is incoherent on all channels.

Relaxing the noise-free assumption:

• (technically should incorporate noise around Equation A13b, and carry the noise through the math above -- but for now, just look at Equation tmp2 above, substituting the signal+noise model):

² Technically, probably should have incorporated a noise model earlier and carried the noise terms through the math.

Expanding the tmp3 above, (FOILing the numerator and denominator)

By applying Assumption 1 we zero-out the terms that are a product of signal and noise:

And by applying Assumptions 2 the last term in the numerator also goes to zero (but the last one in the denominator does not!

Leaving us with

But this says that the denominator is biased up by the noise autopower in $H_y$, so $Z_{xy}$ will be smaller than it ought to be.

Other channel combinations will yield biasing up by $H_y$ in the numerator, or by $E$ in the numerator and denominator respectively.

Proof left as an exercise.

Aside: since noise levels in $H$ tend to be lower than in $E$, if you have to have bias it maybe better to have it from H. That said, you can also compute the biased estimates and average them -- this was the recommended approach until the Remote Reference method.

Remote Reference¶

The above approach motivates the RR approach as an unbiased estimator. WLOG we can replace channels $A$ and $B$ in the matrix $C$ with the remote $H_x$ and remote $H_y$ respectivley. If we assume that the remote reference channels are parallel to the channels at our station then all the assumptions 0,1,2 can be applied and we simply wind up with tmp4d, but in this case:

And by Assumption 2 the second term in the denominator is 0 and we recover an unbiased estimate:

This “cross-power trick” is a common enough in geophysics (and other fields) that it is worth looking at -- If we had started from first principles on a noise-free homogeneous earth we could have derived that

• Proof left as an exercise

but when noise is present in the measurements, we can use the “cross-power trick” to get an unbiased estimate ... this basically selects the part of the signal that is coherent with $R_y$. As long as $R_y$ is relatively quiet, and the noise it does see is incoherent with the local station noise, then we are effectively projecting the local channels onto the coherent signal. It can be thought of as is something like

where $cos({\theta})$ is like the angle between the real answer vector and the answer+noise, although for complex numbers the inner product is not quite as geometrically straightforward to sketch. The same trick can be used in controlled source applications for example to reduce noise in the receiver by taking

Discussion:¶

Before moving on, consider a 3-D earth ... then in 13b_with_R, we would not have had much simplification -- we cannot say in general that $Z_{xx}$ is 0 for one thing, so we would never have been able to strike the $H_x$ autopowers like we did in tmp2 -- because we would not have been able to “zero-out” the first terms in the sums of numerator and denominator as was done in tmp0, tmp1.

We would have a bigger, messier algebraic expression to track.

But what about Assumptions 1 & 2.

• Assumptions 1 is a topic for discussion. It would seem that if you have noise that is effectively plane wave, as long as it is present at both local and remote, you should be able to treat is as signal -- it walks like a duck -- (CSMT, CSAMT).
• Assumptions 2 is not the same -- incoherent channel noise is rather optimistic. A ground sqirrel chewing a cable is one wat this can happen, so is a faulty electronic component, but for a real physical process that generates EM Noise, assuming that it effects only one channel is not realistic. For one thing, E, H are coupled, so it is only plausible that noise in an H would have a corresponding noise in E (Discuss electric fences), and also, it woudl be very lucky if a noise source was so precisely polarized that it coupled into only one of say $H_x$ or $H_y$, more likely it would be in both -- and ditto for electrics.

RR Summary:¶

The RR method will still work well for 2 or 3D, BUT STILL ASSUMES that channel noise is uncorrelated with other channel noise to get an unbiased estimate. So-- N.B. If the noise is not incoherent, i.e. say some channels at either site have noise correlated with noise at the same site, we do not get the friendly cancelations of the previous section, rather the noise contaminates the estimates. The details of working through that are somewhat algebraically involved, and it is probably better to build a toy model with coherent noise in various channel combinations and look at the bias rather than derive it -- at least for now.

• For incoherent, Gaussian noise, we get an unbiased estimate of the TF
  • these converge to the true value as the number of observations increases

Practical Notes:¶

• The time series described in Equation 3a are a sequence of observations.
  • Without multiple observations, the solution to Equation 1 is underdetermined.
• With noise in the measurements, incoherent or otherwise, we need to apply statistical methods to get TF estimates
  • this involved compensation for outliers (leverage points)
  • estimates improve as the number of observations increases.
  • Noise sources may be transient and entire blocks of time series may be better discarded than processed
• Processing should also consider nothc filtering power-lines or other nearby sources.

References:¶

Sims, W.E., F.X. Bostick, and H.W. SMITH. “THE ESTIMATION OF MAGNETOTELLURIC IMPEDANCE TENSOR ELEMENTS FROM MEASURED DATA?” GEOPHYSICS 36, no. 5 (1971). Sims et al. (1971).

Vozoff

APPENDIX¶

Notes on Vozoff 1991:¶

Note that $Z_{xx}$ expression 13a is exactly (term-for-term) that of Vozoff Equation 42a.

Note that $Z_{xy}$ is 13b is nearly Vozoff 42b, but some terms are rearranged.

FWIW, The steps to turn 13b into Vozoff42b are: 1: Multiply top and bottom by -1

Then switch order of addition so that postive terms are first

And then switch permute the multiplcations in the denominator

There is one more difference with Vozoff91, which is that V91 Equation 42b has a typo: The very last term in bra-kets in the denominator (above $H_x A^*$) is unfortunately, incorrectly transcribed in V91 as (above $H_y A^*$).

Notes on invertability , and independence of $C$¶

Recall the determinant

so we dont want

i.e. the following is NOT OK

For example, we cannot allow A equal B since if we set A=B

Which is to say the matrix $C$ must have columns that are linear independent.