The command centres the log measurements and performs principal component analysis of their covariance matrix. With cor=F, it does not first standardize the variables to variance one. For a unit principal component loading vector $q$, the score $q^T(z_i-\bar{z})$ has sample variance $q^T\hat{\Sigma }q$. By the spectral theorem for real symmetric matrices, expand $q$ in an orthonormal eigenbasis; that variance is a weighted average of the eigenvalues. Its maximum is the largest eigenvalue, reached at the corresponding eigenvector. Repeating subject to orthogonality to previous directions gives the other principal components.

The software uses the divisor-$31$ sample covariance matrix for its component standard deviations, whereas var used divisor $30$. Diagonalizing the displayed matrix gives approximatelyTaking square roots of the latter gives $(0.5671603,0.07312405,0.03064834)$, agreeing with the output to the rounding of the displayed sample covariance matrix. The normalization changes these numerical standard deviations but leaves the eigenvectors and proportions of explained variance of a principal component unchanged.

The first principal component explains about $98.08\%$ of the total log-variable variance, the second another $1.63\%$, and the last about $0.286\%$. Keeping two preserves about $99.714\%$. An original scree plot displays the sharp drop after the first component:

Consequently the centered cloud is very elongated along one direction: a one-dimensional approximation represents almost all of its squared variation. From the displayed covariance matrix, a unit first eigenvector, with its arbitrary sign chosen positive, is approximatelyIts principal component score is therefore $0.3980(\log G_i-\overline{\log G})+0.09514(\log H_i-\overline{\log H})+0.9124(\log V_i-\overline{\log V})$. All coefficients have the same sign, so this direction describes overall tree size, with log volume dominating the unstandardized variation. The summary alone supplies proportions, not the loading directions; the displayed sample covariance matrix supplies the latter. This is a geometric approximation, not evidence of an exact deterministic relation or independent original variables.

One further use is principal component analysis on a correlation matrix, obtained with cor=T, which equalizes marginal variances and answers a different scaling-sensitive question. Another is principal component regression: fit a response to selected principal component scores rather than highly correlated original predictors, choosing the number using prediction assessment. The fitted transformation also projects new measurements into the same low-dimensional coordinates, provided the training centering and scaling are retained.

Take logarithms entry by entry before centering. For tree $i$, let $z_i=(\log G_i,\log H_i,\log V_i{)}^T$, and put $\bar{z}=31^{-1}{\sum }_i{z}_i$. The command var(log(b)) computes the sample covariance matrix with the unbiased estimator normalizationThus entry $(j,k)$ is the sum of products of the centered log measurements divided by $30$, not the logarithm of the covariance between the original measurements. Its diagonal entries are the three log-variable sample variances, and the off-diagonal entries are their sample covariances. The positive off-diagonal entries indicate positive association of the log measurements. Taking logarithms turns a multiplicative relation such as $V\approx c{G}^{2}H$ into the approximately linear relation $\log V\approx \log c+2\log G+\log H$, which motivates investigating principal components of these transformed observations.

The output describes complete-linkage clustering: the distance between two current groups is the largest original Jaccard distance between a member of one and a member of the other. Starting with singleton groups, agglomerative hierarchical clustering merges a closest pair at each step. In the merge array, a negative index denotes the original row of that number, and a positive index denotes the group created at that earlier merge. A printed height is the dissimilarity at which the new group forms.

Reading those references recursively gives the following dendrogram construction, with exact heights recovered from the binary profiles:

For example, a cut strictly between $1/2$ and $5/8$ gives four groups: Chad–Nicolas–Frederic–John–Garfield–Fred–Gbenga, Graham–Tim–Philip–Mark, Sauli, and Juliet. A cut strictly between $3/4$ and $7/8$ gives two larger groups. The leaf order is not intrinsic; the merge memberships and heights are. Reflecting branches leaves the dendrogram unchanged in meaning. The zero-height pair indicates identical recorded responses, while high joins indicate comparatively dissimilar profiles. These conclusions depend on the choice of Jaccard distance, which ignores shared negative answers; the dendrogram itself does not establish substantive demographic classes.

Treat each row as a set $A_i$ of attributes with response one. The binary option computes the Jaccard distance, rather than counting matches at joint absences. For a pair of rows, let $a$ count coordinates equal to one in both, $b$ count $(1,0)$, $c$ count $(0,1)$, and $d$ count $(0,0)$. ThenThe denominator counts coordinates present in at least one row; joint absences do not contribute. If both rows have no presences, use distance zero. This differs from normalized Hamming distance, which would divide $b+c$ by all ten coordinates.

For Philip and Chad, four attributes are present in both and eight in at least one, giving Jaccard distance $4/8=0.50$. Graham and Tim differ in one of eight coordinates present in their union, giving $1/8=0.125$, printed as $0.12$. Fred and Gbenga have identical profiles, hence distance zero despite being different students. The resulting dissimilarity matrix is symmetric with zero diagonal; dist2full restores the other half and the diagonal from the stored pairwise entries. The display rounds entries for readability; the agglomerative hierarchical clustering uses the original unrounded distances.

Under the unrestricted model, insert the maximum-likelihood estimates $(\bar{x},S)$. Under the diagonal-covariance restriction, the quadratic form is minimized over the mean at $\bar{x}$, and the remaining objective separates over the diagonal entries:Differentiating in $v_j>0$ gives ${\hat{v}}_j=S_{jj}$; the derivative changes from negative to positive there. Put $D=\mathrm{diag}(S_{11},\dots ,S_{pp})$. The fitted trace term is $p$ under either model, since $\mathrm{tr}(D^{-1}S)=p=\mathrm{tr}(S^{-1}S)$. Consequently the likelihood-ratio test statistic is determined byThis is the Gaussian diagonal-covariance likelihood-ratio test. Reject for sufficiently small $\log |R|$, equivalently for large $-n\log |R|$. Under the regular large-sample calibration, Wilks theorem gives $-n\log |R|\Rightarrow {\chi }_{p(p-1)/2}^2$, since diagonality removes that many free covariance parameters. The exact critical constant can instead be obtained from the null distribution at the actual sample size.

For $p=2$, write the sample correlation coefficient as $r$. Then $|R|=1-r^2$, so the likelihood-ratio test rejects for large $|r|$, not merely for large positive $r$. Under the null, the two variables have a bivariate normal distribution and are independent, and the exact test isIndeed, condition on the centered first-variable observation vector. Decompose the second centered Gaussian vector into its component along that direction and the $n-2$-dimensional orthogonal complement. Their squared lengths give independent ${\chi }_1^2$ and ${\chi }_{n-2}^2$ variables. Their ratio shows that $r\sqrt{(n-2)/(1-r^2)}$ has Student's t-distribution with $n-2$ statistical degrees of freedom.

With the divisor-$n$ sample covariance matrix defined above, the maximum-likelihood estimators areTheir joint sampling distribution is specified byHere $W_p(\nu ,V)$ uses the Wishart distribution convention $\mathbb{E}W=\nu V$. Thus the covariance maximum-likelihood estimator has expected value $(n-1)V/n$; the unbiased estimator of the covariance matrix is $nS/(n-1)$. These distributions follow from orthogonal projections of the jointly multivariate normal observations onto the constant-observation direction and its orthogonal complement, though no proof is required here. The covariance maximum in the positive-definite parameter space requires $S$ nonsingular, which holds almost surely when $n>p$ and $V$ is positive definite. Otherwise the rank of a centered sample covariance matrix shows why the usual unrestricted finite maximum does not exist.

Use the divisor-$n$ sample covariance matrix, since this is the normalization that enters the multivariate normal likelihood function:For a covariance matrix that is a positive-definite matrix $V$, the product of the multivariate normal densities givesTo separate the sample mean from the sample covariance matrix, expand $x_j-\mu =(x_j-\bar{x})+(\bar{x}-\mu )$. The cross terms sum to zero because ${\sum }_j(x_j-\bar{x})=0$. Also $u^T{V}^{-1}u=\mathrm{tr}(V^{-1}u{u}^T)$ by the trace identity. ThereforeThe printed expression suppresses the parameter-independent constant $np\log (2\pi )$; dropping this constant has no effect on maximum likelihood estimation or the likelihood-ratio test.

A curved exponential family restricts a higher-dimensional full canonical family to a lower-dimensional smooth parameter surface:The derivative matrix $J(\theta )=\partial \eta /\partial \theta$ has rank $q$. A genuinely curved surface is not merely an affine subfamily under a choice of natural coordinates. With $n$ iid observations the ambient sufficient statistic is $T_n={\sum }_{i}t(X_i)$, and the parameter score and information areThe second-derivative terms from $\eta (\theta )$ have zero expectation because the centered sufficient statistic has mean zero. Under regularity, the constrained MLE is still first-order asymptotically normal with this information. Curvature affects higher-order bias, likelihood geometry and efficiency; concavity of the ambient likelihood need not survive along the nonlinear parameter surface.

For example, $N(\mu ,{\mu }^2)$ with $\mu >0$ has $t(x)=(x,x^2)$ and natural parameters ${\eta }_1=1/\mu$, ${\eta }_2=-1/(2{\mu }^2)$, constrained by ${\eta }_2=-{\eta }_1^2/2$. It has one free parameter but generally still needs both $\sum X_i$ and $\sum X_i^2$ as sufficient statistics. To see why, a likelihood ratio between two samples is constant in $\mu$ only if their statistic differences satisfy $\Delta T_1/\mu -\Delta T_2/(2{\mu }^2)=\text{constant}$ for all positive $\mu$, forcing both differences to zero. Lower parameter dimension therefore does not generally reduce the minimal sufficient statistic to that dimension.

Several finite-sample counterparts of the influence function should be distinguished. The empirical influence function substitutes the empirical distribution: ${\mathrm{EIF}}_n(x)=\mathrm{IF}(x;T,F_n)$. For a smooth M-estimator it isIt requires a well-defined derivative and nonsingular denominator. An unsmoothed empirical quantile functional, for example, need not admit the population-style density derivative.

The actual one-observation sensitivity curve is a finite difference:It records contamination of size $1/(n+1)$ rather than an infinitesimal perturbation. For a smooth functional it approaches the population influence curve as the sample grows. A deletion counterpart uses $F_{n,-i}=(n{F}_n-{\delta }_{X_i})/(n-1)$, givingThe jackknife resampling pseudovalue $J_i=nT(F_n)-(n-1)T(F_{n,-i})$ is therefore approximately $T(F_n)+{\mathrm{EIF}}_n(X_i)$. The jackknife variance estimate $(n-1)n^{-1}{\sum }_i[T(F_{n,-i})-{\overline{T}}_{(-)}{]}^2$ corresponds asymptotically to $n^{-2}{\sum }_i{\mathrm{EIF}}_n(X_i{)}^2$ for centered empirical influences.

In linear regression, Cook's distance measures the fitted-value effect of deleting a case. For $p$ coefficients, residual $e_i$, leverage $h_{ii}$ and variance estimate ${\hat{\sigma }}^2$, its usual expression is $D_i=e_i^2{h}_{ii}/[p{\hat{\sigma }}^2(1-h_{ii}{)}^2]$. It illustrates that a case can be influential through both a large residual and high leverage. Finite perturbations and nonlinear refitting can differ materially from a first-order influence approximation; these diagnostics should not be conflated.

A one-sample U-statistic of fixed degree $m$ has a symmetric U-statistic kernel $h$ and averages it over all distinct subsets of an iid sample:It is unbiased for $\theta =\mathbb{E}h(X_1,\dots ,X_m)$. With finite second moment, define the first Hoeffding projection $h_1(x)=\mathbb{E}[h(x,X_2,\dots ,X_m)]-\theta$. Orthogonality of the higher projections givesIf ${\zeta }_1=\mathrm{Var}(h_1(X))>0$, the central limit theorem therefore proves the U-statistic central limit theoremEquivalently the influence function of a U-statistic is $m{h}_1(x)$, obtained by differentiating its population functional $\int hd{F}^m$ in each distribution argument. If ${\zeta }_1=0$, the root-$n$ limit degenerates and the second or higher projection determines a different, often nonnormal, limit.

For example, $h(x,y)=(x-y{)}^2/2$ gives the unbiased sample variance with denominator $n-1$. Its first projection is $[(x-\mu {)}^2-{\sigma }^2]/2$, so the influence function is $(x-\mu {)}^2-{\sigma }^2$. When the fourth moment is finite, the leading variance of the sample variance is $(\mathbb{E}(X-\mu {)}^4-{\sigma }^4)/n$. This example connects unbiased symmetric estimation to the functional influence calculation.

Partition a regular parameter into interest $\psi$ and nuisance $\lambda$. Parameter orthogonality means the cross-block of the expected Fisher information matrix vanishes:It is expected-information orthogonality, not a demand that every observed cross derivative be zero. To first order, the asymptotic normal distribution of the joint maximum-likelihood estimator has independent interest and nuisance blocks. Its efficient interest information simplifies from the Schur complement $I_{\psi \psi }-I_{\psi \lambda }{I}_{\lambda \lambda }^{-1}{I}_{\lambda \psi }$ to $I_{\psi \psi }$.

For scalar interest, new nuisance coordinates $\eta$ can locally be chosen by the local orthogonal nuisance reparametrizationIndeed the new interest score is $U_{\psi }+({\partial }_{\psi }\lambda {)}^{\top }{U}_{\lambda }$, and its expected product with the transformed nuisance score is zero by this equation. Smoothness, nonsingular nuisance information and suitable initial coordinates provide the local ordinary-differential-equation construction. Simultaneous orthogonalization for several interest coordinates may face integrability constraints.

For normal location and variance $(\mu ,{\sigma }^2)$, the cross information is zero: the location score is proportional to the centered observation and the variance score to its squared value minus its expectation, whose product has zero expectation by normal symmetry. This is a useful orthogonal parameterization. In general orthogonality reduces nuisance sensitivity and motivates adjusted profile constructions, but it implies neither exact estimator independence nor an exact nuisance-free conditional likelihood.

In a transformation model, a group acts on data $y$ and parameters compatibly: $Y\sim P_{\theta }$ implies $gY\sim P_{g\theta }$. A statistic $A$ is invariant if $A(gy)=A(y)$, and is a maximal invariant if equality of its values also implies the samples lie on the same group orbit. A maximal invariant records all information unchanged by the transformations; measurable invariant procedures therefore factor through it, subject to the ordinary orbit-space measurability conditions.

An equivariant estimator satisfies $T(gy)=gT(y)$ for the induced action on the quantity estimated. For the positive affine location-scale group, $gy=a+by$ with $b>0$, this means $T_{\mathrm{loc}}(a+by)=a+b{T}_{\mathrm{loc}}(y)$ and $T_{\mathrm{scale}}(a+by)=b{T}_{\mathrm{scale}}(y)$.

For a nonconstant ordered sample, take $\bar{y}$ and $s^2=n^{-1}{\sum }_i(y_i-\bar{y}{)}^2>0$. Thenis a maximal invariant: it is unchanged by the positive affine group, and if $A(y)=A(z)$ then $z_i=\bar{z}+(s_z/s_y)(y_i-\bar{y})$, explicitly recovering the transformation. Both mean and standard deviation are equivariant. If the whole model varies only by location and positive scale, the normalized sample's distribution is independent of those two parameters, giving an ancillary statistic. It may still retain information about additional shape parameters.

Equivariance yields a coherent transformation rule and, under an invariant loss, can simplify risk comparison to a reference parameter. It does not on its own establish optimality. Normalization by a parameter-space estimator with a nontrivial stabilizer requires quotienting that remaining action; otherwise an apparently standardized statistic need not be maximal.

The Laplace approximation evaluates an integral dominated by an interior maximum. Suppose $h$ has a unique interior maximizer $x_0\in {\mathbb{R}}^d$, negative-definite Hessian there, and sufficient smoothness and tail control. Put $H=-{\nabla }^{2}h(x_0)>0$. Taylor expansion under $x=x_0+u/\sqrt{n}$ gives a Gaussian leading integrand andThe leading factor follows from the multivariate Gaussian integral. Odd cubic terms integrate to zero over the limiting symmetric domain; fourth-order terms and squared cubic terms contribute to the first correction, assuming enough derivatives. If $g(x_0)=0$, a different leading term may be needed.

For a likelihood integral with prior density $\pi$, a common version is $\int e^{\ell (\theta )}\pi (\theta )d\theta \approx e^{\ell (\hat{\theta })}\pi (\hat{\theta })(2\pi {)}^{d/2}|j(\hat{\theta }){|}^{-1/2}$, because observed information is of order $n$. This is useful for posterior normalization, marginal likelihood and approximating integrated nuisance effects. A prior or integration measure is part of that construction; integrating nuisance parameters is different from merely profiling them. Boundary maxima, degenerate curvature, multiple modes or uncontrolled tails invalidate the simple one-mode formula.

For parameter of interest $\psi$ and nuisance parameter $\lambda$, the profile likelihood has log form ${\ell }_p(\psi )=\ell (\psi ,{\hat{\lambda }}_{\psi })$, where ${\hat{\lambda }}_{\psi }$ is the constrained nuisance maximum. Profiling optimizes over the nuisance dimension but does not account for its local volume or conditioning geometry. The modified profile likelihood corrects this loss.

Use ancillary sample coordinates $(\hat{\psi },\hat{\lambda },a)$, and let $j_{\lambda \lambda }(\psi ,{\hat{\lambda }}_{\psi })=-{\ell }_{\lambda \lambda }$ be the nuisance block of observed information. Holding $\hat{\psi },a$ fixed, put $D(\psi )=|\partial {\hat{\lambda }}_{\psi }/\partial \hat{\lambda }|$. The Barndorff-Nielsen modified profile likelihood is, up to data-only factors,This combines a nuisance-information determinant and a data-coordinate Jacobian determinant. Implicit differentiation of ${\ell }_{\lambda }(\psi ,{\hat{\lambda }}_{\psi })=0$ gives $\partial {\hat{\lambda }}_{\psi }/\partial \hat{\lambda }=j_{\lambda \lambda }^{-1}{\ell }_{\lambda ;\hat{\lambda }}$, so equivalentlyThe semicolon denotes differentiation in the varying sample coordinate, rather than a parameter derivative. This form also explains the positive adjustment in Question 2: in that canonical conditional problem the mixed sample derivative is data-only, leaving the positive half-log curvature. The adjustment sign cannot be copied between unrelated constructions without tracking their Jacobians.

In an appropriate orthogonal nuisance formulation, the Cox-Reid adjusted profile likelihood uses ${\ell }_p-\frac{1}{2}\log |j_{\lambda \lambda }|$. It is a useful simpler adjustment, but discarding $D$ is not universally legitimate. For normal data with interest ${\sigma }^2$ and nuisance $\mu$, the constrained and full nuisance fits both equal $\bar{X}$, so $D=1$ and $j_{\mu \mu }=n/{\sigma }^2$. Thusrecovering the residual degrees of freedom and fitted variance ${\sum }_i(X_i-\bar{X}{)}^2/(n-1)$. These adjustments target nuisance-related bias and higher-order likelihood inference; they are approximations to conditional or marginal constructions, not an assertion that profiling has become an exact full sampling density.

A Bartlett correction improves calibration of a regular likelihood-ratio test statistic $W=2\{\ell (\hat{\theta })-\ell ({\hat{\theta }}_0)\}$, whose null limit is ${\chi }_q^2$ for $q$ tested restrictions. Ifthe corrected statistic is $\boxed{W_B=W/(1+b/n)}$. Its null mean is $q+O(n^{-2})$ and it is compared with the same limiting chi-squared distribution. The coefficient is model-specific, may involve nuisance parameters, and can be calculated from likelihood cumulants or a suitable null simulation; substitution must respect the asymptotic regime.

Under the regular higher-order likelihood expansions and identities that establish Bartlett correctability, this rescaling removes the leading $O(n^{-1})$ distributional error and can give $O(n^{-2})$ chi-square calibration error. Matching one moment by itself does not prove that assertion for an arbitrary statistic, and nonregular boundary or singular models need different treatment.

For illustration, testing a normal variance with unknown mean gives $W=n[-\log (U/n)+U/n-1]$ with $U\sim {\chi }_{n-1}^2$. The gamma logarithmic moment yieldsThus the Bartlett correction for a normal variance with unknown mean uses $b=11/6$, not the correction for the different known-mean model. This example also distinguishes the relative coefficient $b$ from the full divisor $1+b/n$.

For the Wilcoxon signed-rank test, assume independent observations from a continuous distribution symmetric about the specified ${\theta }_0$. Form $D_i=X_i-{\theta }_0$, rank the absolute differences $|D_i|$, and add the ranks of positive differences:Under the continuous null there are no zeros or tied absolute values with positive probability. Symmetry makes the signs independent fair signs conditional on the magnitudes, so the exact null distribution is that of ${\sum }_{j=1}^{n}j{B}_j$ for independent $B_j\sim \mathrm{Bernoulli}(1/2)$. In particularThe coefficients of this polynomial give exact tail probabilities, or they can be computed by recursively adding a rank with probability one half. A two-sided test rejects for unusually small $\min (W_{+},W_{-})$, equivalently for extreme $W_{+}$ on either side of its center. Choose symmetric tail cutoffs with total null probability at most $\alpha$; discrete boundary randomization can attain size exactly $\alpha$. An exact conservative two-sided p-value is $\min \{1,2\min ({\mathbb{P}}_0(W_{+}\le w),{\mathbb{P}}_0(W_{+}\ge w))\}$. For a positive-location alternative use the upper tail, and for a negative-location alternative use the lower tail.

For large samples, centering and scaling by the displayed mean and standard deviation gives an asymptotic standard normal statistic. This follows from the independent weighted Bernoulli representation, because the largest squared weight is negligible relative to the sum of squared weights. A continuity correction improves the normal tail approximation. With observed zero differences, omit them in the usual convention; tied absolute differences receive average ranks, and the conditional null variance is then one quarter of the sum of squared assigned ranks, not automatically the no-tie formula.

Symmetry is essential: having median ${\theta }_0$ alone need not make the signs independent of the magnitudes, and so does not justify these exact critical values. The test is particularly useful for location shifts in symmetric families. It is a valid test of the stated symmetry null, but is not an omnibus procedure guaranteed to detect every asymmetric alternative.

For $n\ge 2$, the probability integral transform makes $U_i=F(X_i)$ iid uniform on $[0,1]$. Continuity of $F$ is enough; strict increase is unnecessary. Their extreme order statistics have joint densitychoose the observations giving the two extremes and put the remaining $n-2$ between them. Transform to $r=v-u$ and integrate $u$ from zero to $1-r$. The uniformized sample range therefore hasEquivalently $R\sim \mathrm{Beta}(n-1,2)$. For $n=1$, the range is identically zero instead.

The exact failure probability isThus the exact integer sample size is the least $n\ge 2$ making this quantity at most $\delta$. For $\epsilon \downarrow 0$ with $n\epsilon \to \vartheta >0$, it tends to $(1+\vartheta )e^{-\vartheta }$. This function decreases strictly from one to zero on positive $\vartheta$, because its derivative is $-\vartheta e^{-\vartheta }$. For $0<\delta <1$ there is therefore exactly one positive solution ofThis gives the extreme-order-statistic tolerance interval sample-size ruleWhen $\delta$ is small, $\vartheta =\log (1/\delta )+\log (1+\vartheta )$, so a useful first refinement is $\vartheta \approx \log (1/\delta )+\log \log (1/\delta )$. This is an approximation, not an exact integer guarantee; use the preceding exact failure probability to check rounding and very stringent tail requirements. The interval here encloses population content, rather than estimating a fixed distribution parameter.