Skip to content

Home / Statistical Tools / MSA / Extended / Math Details

Math Details

The mathematics behind an Extended measurement system analysis. An Extended study uses the expected mean squares method or restricted maximum likelihood; the XbarR method is not available, so it is not covered here.

Everything on this page applies to the Crossed and Nested studies too. What is particular to Extended is that the model is declared rather than fixed, which brings in fixed terms, deeper nesting, and the model that restricted maximum likelihood cannot fit.

The model and the components

A measurement system analysis fits a variance components model. Every measured value is written as the sum of the study's overall mean, one contribution per model term, and a residual:

\[ y = \mu + \sum_{t} \tau_t + \varepsilon \]

\(\mu\) is the grand mean, \(\tau_t\) is the contribution of model term \(t\) (a part effect, an operator effect, an interaction, or an additional factor in an Extended study), and \(\varepsilon\) is the residual. Each random term is taken to vary with its own variance component, and the residual with \(\sigma^2_{\text{repeat}}\):

\[ \operatorname{Var}(\tau_t) = \sigma^2_t, \qquad \operatorname{Var}(\varepsilon) = \sigma^2_{\text{repeat}} \]

Estimating those \(\sigma^2_t\) is the whole job of the analysis. Which estimation method is used decides how they are computed, but not what they mean, and not how they are then grouped and reported.

In an Extended study the terms are whatever the Model tab declares: one main effect per factor, each interaction left checked, and the residual. Nesting changes what a term's levels are rather than adding a term.

Every estimated component is placed into exactly one of three buckets. The rule has three clauses and no exceptions:

Bucket Which components go in it
Repeatability the residual, and nothing else
Part to part every term selected under Part-to-part variation
Reproducibility every other term in the model

So reproducibility is defined by what is left over rather than by a list of its own. An interaction, an additional factor in an Extended study, and a term whose estimate came from a fixed effect all roll into reproducibility unless the user moved them to part to part.

That last point is what the Part-to-part variation checkbox list on the Model tab does, and the dialog says so in the same words: by default the variation from additional factors and interactions is attributed to reproducibility, and selecting a term attributes its variation to part to part instead.

The Part-to-part variation group on the Model tab is what moves a term between the second and third buckets of that rule, and it is the only place a term's bucket can be changed.

The three buckets are summed, and two further aggregates are formed from those sums:

\[ \sigma^2_{\text{repeat}} = \text{the residual component} \]
\[ \sigma^2_{\text{reprod}} = \sum_{t \,\notin\, \text{part}} \sigma^2_t \]
\[ \sigma^2_{\text{part}} = \sum_{t \,\in\, \text{part}} \sigma^2_t \]
\[ \sigma^2_{\text{gage}} = \sigma^2_{\text{repeat}} + \sigma^2_{\text{reprod}} \]
\[ \sigma^2_{\text{total}} = \sigma^2_{\text{gage}} + \sigma^2_{\text{part}} \]

\(\sigma^2_{\text{gage}}\) is reported as Total Gage R&R and \(\sigma^2_{\text{total}}\) as Total Variation. Every percentage on the report is a ratio against one of these two, so they are computed once and never recomputed further down the sheet.

The report prints one row per bucket total and one row per member component, in this fixed order: Total Gage R&R, Repeatability, Reproducibility total, each reproducibility member, Part to part total, each part to part member, Total Variation. Members appear in ascending model term order.

Full details: Variance Components.

Fixed terms

A factor declared Fixed on the Model tab has no variance component. What the report shows in its place is defined here.

A term whose factors are all declared Fixed has no variance component, because a fixed effect has no distribution to have a variance. So that it can still take a place on the report, a substitute quantity is reported in the component column: the sum of the term's squared coefficients divided by its number of levels.

\[ \text{substitute}_t = \frac{1}{J_t}\sum_{j=1}^{J_t} \alpha_{tj}^2 \]

\(\alpha_{tj}\) are the term's estimated coefficients and \(J_t\) its number of levels. The row is marked as a fixed term substitute, and two consequences follow that a reader will notice:

  • It is not a variance component, even though it sits in the variance component column and is summed into an aggregate along with the real ones.
  • It carries no confidence interval, because a fixed effect has no sampling distribution for an interval to be built from.

What is reported for each component

Each reported component is also reported as a standard deviation, which is simply its square root:

\[ \sigma_t = \sqrt{\sigma^2_t} \]

Because the floor at zero is applied to the variance before the square root is taken, this is always defined.

The study variation of a component is its standard deviation multiplied by \(k\), the study variation multiplier:

\[ SV_t = k\,\sigma_t \]

\(k\) is set by Study variation multiplier (k): on the Extended options group and defaults to 6, which spans a range of plus and minus three standard deviations. A value of zero or below selects the default rather than being used.

Percent contribution compares variances:

\[ \%\text{Contribution}_t = 100 \times \frac{\sigma^2_t}{\sigma^2_{\text{total}}} \]

Because the parts are variances and variances add, these percentages sum to 100 across the components.

Percent study variation compares standard deviations:

\[ \%\text{Study Variation}_t = 100 \times \frac{\sigma_t}{\sigma_{\text{total}}} \]

The multiplier \(k\) cancels, so this is the same number whether it is formed from standard deviations or from study variations.

On the report this quantity is not a column. It appears as four named rows beneath the results table, one for each aggregate, under the abbreviations the gage literature uses:

Row Which component's percent study variation
%GRR Total Gage R&R
%EV Repeatability, the equipment variation
%AV Reproducibility, the appraiser variation
%PV Part to part, the part variation

A row whose component the study did not produce reports NA. The member components of reproducibility and of part to part have no row of their own here, so their percent study variation is not reported at all even though the aggregate ones are.

These two percentage columns answer different questions and do not agree, and neither is wrong. Percent contribution is a ratio of variances; percent study variation is a ratio of standard deviations, which is the square root of the same ratio. A component contributing 9 percent of the variance contributes 30 percent of the standard deviation. Percent study variation does not sum to 100 across the components, because standard deviations do not add.

Percent of tolerance compares a component's study variation against the tolerance width:

\[ \%\text{Tolerance}_t = 100 \times \frac{SV_t}{\text{Tolerance}} = 100 \times \frac{k\,\sigma_t}{\text{Tolerance}} \]

Unlike the two percentage columns above, this one does move with \(k\). It is not reported when no specification limit was supplied.

When a historical standard deviation is supplied, each component is also reported against it:

\[ \%\text{Historical Process}_t = 100 \times \frac{\sigma_t}{\sigma_{\text{hist}}} \]

This column is not reported when no historical value was supplied. Note that it uses the historical value directly as the denominator, with no subtraction; that is different from how the historical value enters the misclassification calculation.

The precision to tolerance ratio compares the whole measurement system against the tolerance:

\[ P/T = \frac{k\,\sigma_{\text{gage}}}{\text{Tolerance}} \]

It is not reported when no specification limit was supplied.

The number of distinct categories estimates how many separate groups of parts the measurement system can tell apart:

\[ ndc = \left\lfloor 1.41 \times \frac{\sigma_{\text{part}}}{\sigma_{\text{gage}}} \right\rfloor, \qquad ndc \geq 1 \]

The ratio is truncated toward zero and then floored at 1, so the reported count is never zero and never negative. A ratio large enough to exceed the integer range is clamped to the largest reportable integer.

Two properties of this number are worth stating because they are easy to assume otherwise:

  • There is exactly one of it. It is always computed from the estimated part to part and total gage standard deviations, even when a historical standard deviation was supplied. No historical variant of the count exists.
  • It is a truncation, so a reported count of 4 covers a range of underlying ratios. Where a confidence interval on the count is reported, that interval is on the underlying continuous ratio, not on the truncated integer.

The count is not reported at all when the total gage standard deviation is exactly zero, because the ratio it is a truncation of would divide by zero. The report prints a sentence in place of the number rather than a blank cell or a zero, since a zero would read as a measured count of zero categories.

The test is on total gage R&R, not on repeatability alone. A study in which every repeated measurement repeats its own value exactly has no repeatability variation, but if its operators still disagree it has a real gage variation, and its count is reported.

Full details: Gage Evaluation.

How the components are estimated

Expected mean squares

The expected mean squares method builds the analysis of variance table, works out what each mean square is expected to be in terms of the variance components, and then solves that system for the components.

Its name is what it does. Every mean square in the table estimates some combination of the variance components; writing down those combinations gives a set of equations, and the components are what solves them.

On balanced data this is the classical analysis of variance estimate, which is what the option label means by "Equivalent to ANOVA when balanced".

For each term \(T\) in the model, the mean square \(MS_T\) has an expected value that is a weighted sum of the variance components plus the residual variance:

\[ E[MS_T] = \sigma^2_{\text{repeat}} + \sum_{B} c_{TB}\,\sigma^2_B \]

\(B\) runs over the random terms of the model and \(c_{TB}\) is the coefficient of term \(B\)'s component in term \(T\)'s expected mean square. Collecting these for every term gives a triangular system of equations in the unknown components.

A fixed term's column is left empty. A fixed effect has no variance component, so no coefficient is computed for it and it takes no part in the solve.

When the data are balanced the coefficients have a closed form. The coefficient of \(B\) in \(T\)'s expected mean square is the total number of observations divided by the product of the level counts of the factors in \(B\):

\[ c_{TB} = \frac{N}{\prod_{f \in B} L_f} \]

\(N\) is the number of observations used and \(L_f\) the number of levels of factor \(f\). The coefficient is non-zero only when term \(B\) contains term \(T\).

Read plainly: a component's coefficient is the number of observations that share each of its levels. In a crossed study with \(p\) parts, \(o\) operators and \(r\) replicates, the part component appears in the part mean square with coefficient \(o\,r\), because each part was measured \(o\,r\) times.

The closed form is used only when the data are genuinely balanced, and that is checked rather than assumed.

When the data are not balanced the closed form is wrong, and the coefficients are computed from the model matrices instead. For term \(T\) the projector onto that term's own space is formed,

\[ \mathbf{P}_T = \mathbf{X}\mathbf{G}\mathbf{L}^{\top}\left(\mathbf{L}\mathbf{G}\mathbf{L}^{\top}\right)^{-1}\mathbf{L}\mathbf{G}\mathbf{X}^{\top} \]

where \(\mathbf{X}\) is the model matrix, \(\mathbf{G}\) a generalized inverse of \(\mathbf{X}^{\top}\mathbf{X}\), and \(\mathbf{L}\) selects the columns belonging to \(T\). Each coefficient is then a trace against the other term's own incidence matrix, divided by the term's degrees of freedom:

\[ c_{TB} = \frac{\operatorname{trace}\!\left(\mathbf{P}_T \mathbf{Z}_B \mathbf{Z}_B^{\top}\right)}{\nu_T} \]

\(\mathbf{Z}_B\) is the zero-or-one incidence matrix of term \(B\) and \(\nu_T\) the degrees of freedom of term \(T\).

This is the same construction whatever the design: crossed, nested, or any Extended shape. There is one algorithm, and the balanced closed form above is a shortcut into it rather than a separate method.

The construction carries its own check. If the projector's trace disagrees with the degrees of freedom the analysis of variance reports for that same term, the analysis is refused rather than reporting coefficients that cannot be trusted.

The system is triangular, so it is solved by back-substitution from the bottom up. Each term's component is solved only after every component its own expected mean square depends on has been solved.

The residual mean square starts it off, since \(E[MS_{\text{residual}}] = \sigma^2_{\text{repeat}}\) with no other term in it, so the residual mean square estimates repeatability directly:

\[ \hat{\sigma}^2_{\text{repeat}} = MS_{\text{residual}} \]

Where the lowest term in the model is itself the lowest estimable level, meaning there is one observation per cell and no residual, there is no residual mean square to start from and the error variance enters the solve as zero instead.

Every component therefore ends up as a linear combination of the mean squares, and the weights of that combination are kept. They are what lets the confidence interval layer recognise whether a component came from a single mean square, from a sum, or from a difference, which decides which interval formula applies.

The subtraction in that formula is why this method can produce a negative variance estimate. If the mean square being tested comes out smaller than its own denominator, the difference is negative, and dividing by a positive coefficient keeps it negative.

Nothing about the arithmetic is wrong when that happens. A variance component is a non-negative quantity being estimated by a difference of two noisy numbers, and a difference can land on the wrong side of zero when the true component is small.

The report handles it by reporting zero, keeping the unfloored estimate beside it, flagging the row, and printing this footnote:

A negative variance component estimate was set to zero. Consider the REML method, which estimates variance components under the constraint that they cannot be negative.

Full details: Expected Mean Squares.

Restricted maximum likelihood

Restricted maximum likelihood, abbreviated REML, estimates every variance component in the model at once, from the data as a whole, rather than solving for them from an analysis of variance table. It is the method the Automatic setting chooses when the data are not balanced.

The terms that enter the fit are exactly the terms of the declared model, and nothing is added or dropped by the method itself:

  • one term per declared factor main effect
  • each interaction left checked on the Model tab
  • each nesting relationship declared on the Model tab, which changes what a term's levels are rather than adding a term
  • the residual

A factor declared Fixed enters as a fixed effect, so it contributes to \(\mathbf{X}\) and has no variance component estimated for it. At least one term must be random. A model in which every factor is declared fixed leaves REML with no variance to estimate, and the engine refuses it with a named reason rather than returning zeros; the expected mean squares method reports the same study.

The components are estimated under the constraint that a variance cannot be negative. This is the substantive difference between REML and the expected mean squares method, and it shows up directly on the report:

  • Under expected mean squares a component's moment estimate can come out negative, and zero is reported in its place with a footnote saying so.
  • Under REML nothing negative is ever computed. A reported zero means the constrained estimate reached the zero boundary, and the footnote says that instead.

The product itself points a reader from the first case to the second: the expected mean squares footnote suggests considering the REML method precisely because it constrains the components to be non-negative.

REML cannot fit a model containing a term that involves more than three factors once nesting is counted. Nesting counts here because a nested term carries its parents with it, so a term that looks like two factors in the dialog can involve four once its nesting chain is followed.

When the declared model crosses that limit, the analysis uses the expected mean squares method instead, whether REML was chosen directly or reached through Automatic. The Model tab warns about this before the analysis is run, in these words:

A term in this model involves more than three factors once nesting is counted, so the REML method cannot fit it; if REML is selected or chosen automatically, the analysis will use expected mean squares instead.

The Estimation method: row of the report states the method that was actually used, so the substitution is visible afterwards as well as being warned about beforehand.

Full details, including what is maximised and what is reported about the fit: Restricted Maximum Likelihood.

The F tests and the removal loop

Each term is tested with an F statistic, and the denominator is chosen the same way the component solve chooses its subtraction: the term whose expected mean square equals the numerator's with the tested component removed.

\[ F_T = \frac{MS_T}{MS_{\text{denominator}}} \]

Where that reduced expectation is the residual variance alone, the residual mean square is the denominator. This is an exact F test, and the analysis of variance table says so by not flagging it.

In a more complicated model no single mean square matches what the denominator needs to be. A denominator is then synthesized as a combination of several mean squares whose expected value is the right one, and its degrees of freedom are approximated:

\[ \nu_{\text{denominator}} = \frac{\left(\sum_i a_i MS_i\right)^2}{\sum_i \dfrac{\left(a_i MS_i\right)^2}{\nu_i}} \]

The weights \(a_i\) come from the same elimination that chooses the denominator, and the residual always carries a weight of one minus the sum of the term weights, because every row of the system has an error coefficient of exactly one. A mean square whose weight works out to zero is left out of the printed synthesis rather than listed with a zero beside it, and the residual appears in that list under the name Error.

A negative weight is not a problem for the degrees of freedom: each contribution is squared, so a denominator formed from a difference of mean squares still gets a positive approximation.

Two visible consequences, both of which the table reports rather than hides:

  • The degrees of freedom are not in general a whole number. A fractional denominator degrees of freedom on the report is this approximation, not a rounding error.
  • The test is marked as not exact. The table prints Not an exact F-test. for such a term, and reports the denominator mean square and its degrees of freedom so the test can be traced.

Where the synthesized denominator does not come out strictly positive there is no usable F statistic at all, and the table reports no F and no p for that term rather than an arbitrary number.

When Remove insignificant interactions is on, terms are tested and removed one at a time, and each removed term's degrees of freedom and sum of squares are pooled into the residual. A removed term's variation is therefore counted as repeatability rather than being dropped from the study.

The pooling changes the residual mean square, so every remaining term's F test and every component estimate change with it. Terms are therefore removed in successive passes rather than all at once, and the report records which pass removed a given term.

The report carries two analysis of variance tables when anything was removed, headed ANOVA (All Terms) and ANOVA (Terms used for the Gage R&R calculation), and a removed term is marked (removed) in the first. With nothing removed there is one table.

Probabilities of misclassification

The misclassification block asks what the measurement system does to accept and reject decisions. A part has a true value and a measured value, and they differ by measurement error, so a good part can be rejected and a bad part accepted.

Both quantities are treated as normal. The true value has the study mean and the part to part standard deviation; the measured value differs from the true value by measurement error with the total gage standard deviation:

\[ X \sim N\!\left(\bar{x},\; \sigma^2_{\text{part}}\right), \qquad Y \mid X \sim N\!\left(X,\; \sigma^2_{\text{gage}}\right) \]

\(X\) is the true value and \(Y\) the measured value. "In specification" means lying between whichever specification limits were supplied.

Every probability below is a probability over that pair, so it is a statement about the population of parts and the gage together, not about any part that was actually measured.

The two joint probabilities each cover one way of getting the decision wrong, as a fraction of all parts:

\[ P(\text{good part rejected}) = P\big(X \in \text{spec} \;\text{ and }\; Y \notin \text{spec}\big) \]
\[ P(\text{bad part accepted}) = P\big(X \notin \text{spec} \;\text{ and }\; Y \in \text{spec}\big) \]

Because they are fractions of all parts, they can be added to each other and compared directly against each other.

The two conditional probabilities ask the same two questions of a restricted population: of the parts that really are good, what fraction get rejected, and of the parts that really are bad, what fraction get accepted.

\[ P(\text{reject} \mid \text{good}) = P\big(Y \notin \text{spec} \;\big|\; X \in \text{spec}\big) \]
\[ P(\text{accept} \mid \text{bad}) = P\big(Y \in \text{spec} \;\big|\; X \notin \text{spec}\big) \]

Each is its joint probability divided by the probability of the condition, so the two families carry the same information on different denominators. The conditional numbers are typically much larger than the joint ones, and the difference is entirely the denominator: a conditional probability discards every part that was not in the group it conditions on.

The fifth reported probability is the denominator the conditionals use:

\[ P(\text{part is good}) = P\big(X \in \text{spec}\big) \]

It describes the parts alone and says nothing about the gage. It is reported because without it the conditional probabilities cannot be turned back into fractions of all parts.

When a historical standard deviation is supplied, it replaces the estimated part to part standard deviation in this block only, and it does so by subtraction rather than substitution:

\[ \sigma_{\text{part}} \;\longleftarrow\; \sqrt{\sigma^2_{\text{hist}} - \sigma^2_{\text{gage}}} \]

The historical value is a process standard deviation, meaning it already contains the measurement error, so the gage variance is removed from it to leave the part to part variation. This is why the dialog requires the historical value to be larger than the total gage standard deviation: otherwise the quantity under the square root is not positive.

Three consequences, each of which a reader can see on the report:

  • The report records which standard deviations were actually used, so the substitution is visible.
  • If the historical value fails the requirement, nothing is substituted and nothing is clamped. The block is refused with a reason naming that failure, and the estimated part to part standard deviation is left alone.
  • A historical value that passes the requirement will produce the block even when the estimated part to part variance is zero, which would otherwise have refused it, because the historical path never reads the estimated part variance at all.

The historical value does not enter the number of distinct categories, and it does not change the percent of tolerance columns. It affects this block and the sweep charts.

A conditional probability needs a population to condition on. When the specification limits sit far enough from the study mean, one of those two populations effectively does not exist, and the corresponding conditional probability is withheld rather than reported as a number computed from almost nothing.

The rule: once the distance from the study mean to the nearest specification limit exceeds 8 part to part standard deviations, the conditional probability for the empty population is not reported, and the report names which case it is:

Situation Withheld Still reported
every part is good the probability a bad part is accepted the other conditional, both joint probabilities, and the probability the part is good
every part is bad the probability a good part is rejected the other conditional, both joint probabilities, and the probability the part is good

The rule is symmetric between the two, and 8 standard deviations is a long way out: the fraction of parts beyond it is about \(6.2 \times 10^{-16}\), roughly one part in \(1.6 \times 10^{15}\).

Full details, including the refusals and the sweep: Probabilities of Misclassification.

Bias and linearity

Bias and linearity are the only part of the analysis that needs to know the true value of a part. Everything else in a measurement system analysis compares measurements against each other, so it can say how much the system varies without knowing whether it is right. These two say whether it is right.

That is why both are reported only when a Reference (Optional) column was supplied. Without one the whole block is absent, and its absence is not a problem to be reported.

Bias is the difference between what the system measures and what the part actually is:

\[ b = \bar{y} - x_{\text{ref}} \]

\(\bar{y}\) is the mean of the measurements in the group being reported and \(x_{\text{ref}}\) the reference value. A positive bias means the system reads high.

Each bias row is tested against zero with a t test:

\[ t = \frac{b}{SE(b)}, \qquad \nu = \text{the row's degrees of freedom} \]

and a two-sided p value is read from the t distribution on \(\nu\) degrees of freedom. A confidence interval is reported on the bias itself, at the confidence level set in the options:

\[ b \pm t_{1 - \alpha/2,\; \nu}\; SE(b) \]

The significance level the p value is judged against is the one set by Significance level:, which defaults to 0.05.

Reading the test the way the product does: a p value at or above the significance level means the bias is statistically zero, which is not the same as the bias being zero. It means the study did not have enough evidence to distinguish it from zero.

Linearity asks whether the bias is the same everywhere in the measuring range, or whether it changes as the parts get bigger. It is a straight-line fit of bias against reference value:

\[ b = \beta_0 + \beta_1 x_{\text{ref}} + \varepsilon \]
  • \(\beta_1\), the slope, is how much the bias changes per unit of reference value. A slope of zero means the bias is the same across the range, whatever its size.
  • \(\beta_0\), the intercept, is the bias the line gives at a reference value of zero.

One row is reported per operator plus one pooled row, so an operator whose bias drifts across the range can be told from one whose bias is merely offset.

The fit uses one point per individual measurement, not one point per part average. A study with ten parts measured three times by each of three operators contributes ninety points to the pooled fit, not ten.

Each linearity row reports the fit and then two quantities derived from the slope:

Reported Meaning
intercept and slope each with its standard error, t statistic and p value
R squared the fraction of the variation in bias the line accounts for
regression standard error the spread of the points about the line
residual standard deviation what sets the half width of the confidence band drawn on the chart
degrees of freedom the point count minus two
process variation used the spread of the process the two quantities below are scaled against
\[ \text{Linearity} = \lvert \beta_1 \rvert \times \text{process variation} \]
\[ \%\text{Linearity} = 100 \times \frac{\text{Linearity}}{\text{process variation}} = 100 \times \lvert \beta_1 \rvert \]

Percent linearity reduces exactly to the absolute slope, because the process variation it is divided by is the same one it was multiplied by. That is worth knowing for two reasons: percent linearity does not depend on the process variation at all, and linearity does. The process variation that was used is reported alongside them for exactly that reason.

The confidence interval on bias is reported in an Extended study, at the level set by Confidence level:. The rule that an Extended study reports no confidence intervals applies to the variance component rows, not to the bias table.

Full details: Bias and Linearity.

See Also

  • MSA Math, all seven shared pages
  • Options, the controls that set every quantity above