Home > 04-Estimation > GHM > Functions-GHM > diagnostics_GHM_constraints_multi_d.m

diagnostics_GHM_constraints_multi_d

PURPOSE ^

% Diagnostics Gauss-Helmert model multidimensional

SYNOPSIS ^

function [vgi,Xv,Rgm,nabla_cg,muv,muv1,U1mqc,U2mc]=diagnostics_GHM_constraints_multi_d(r_U,I,d_I,d_G,Am,Bm,Cov_xx,W_xx,W_gg,vv)

DESCRIPTION ^

% Diagnostics Gauss-Helmert model multidimensional

 referring to constraints
 assuming each group of observations belongs to one group of constraints
      number = I, with d_I observations and d_G constraints

  [vgi,Xv,Rgm,nabla_cg,muv,muv1,U1mqc,U2mc]=
            diagnostics_GHM_constraints_multi_d...
   (r_U,I,d_I,d_G,Am,Bm,Cov_xx,W_xx,W_gg,vv)
 
 r_U    = index set for parameters to be estimated (others are nuisance)
 I      = number of observational groups
 d_I    = dimension of observatioanl groups (must be the same for all)
 Am     = G x U Jacobian
 Bm     = G x N Jacobain
 Cov_xx = U x U CovM of all estimated parameters
 W_xx   = U x U WeightM of all estimated parameters
 Wgg    = 
 vv     = N x 1 residuals

 vgi      = I x dI    matrix of residuals
 Xv       = I x 1     vector standardized test statistics (chi^2)
 Rim      = I * d_I^2 redundancy matrices as row vectors 
 nabla_lv = I x 1     lower bound for detectable outliers (maximum eigenvalue)
 muv      = I x 1     sensitivity factor for all parameters
 muv1     = I x 1     sensitivity factor for parameters specified by r_U
 Um1q     = I * d_I^2 diagonal blocks of \bar U_1 (numerator of muv1^2)
 Um2      = I * d_I^2 diagonal blocks of U2


 Wolfgang Förstner
 wfoerstn@uni-bonn.de
 
 See also GaussHelmertModel

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 %% Diagnostics Gauss-Helmert model multidimensional
0002 %
0003 % referring to constraints
0004 % assuming each group of observations belongs to one group of constraints
0005 %      number = I, with d_I observations and d_G constraints
0006 %
0007 %  [vgi,Xv,Rgm,nabla_cg,muv,muv1,U1mqc,U2mc]=
0008 %            diagnostics_GHM_constraints_multi_d...
0009 %   (r_U,I,d_I,d_G,Am,Bm,Cov_xx,W_xx,W_gg,vv)
0010 %
0011 % r_U    = index set for parameters to be estimated (others are nuisance)
0012 % I      = number of observational groups
0013 % d_I    = dimension of observatioanl groups (must be the same for all)
0014 % Am     = G x U Jacobian
0015 % Bm     = G x N Jacobain
0016 % Cov_xx = U x U CovM of all estimated parameters
0017 % W_xx   = U x U WeightM of all estimated parameters
0018 % Wgg    =
0019 % vv     = N x 1 residuals
0020 %
0021 % vgi      = I x dI    matrix of residuals
0022 % Xv       = I x 1     vector standardized test statistics (chi^2)
0023 % Rim      = I * d_I^2 redundancy matrices as row vectors
0024 % nabla_lv = I x 1     lower bound for detectable outliers (maximum eigenvalue)
0025 % muv      = I x 1     sensitivity factor for all parameters
0026 % muv1     = I x 1     sensitivity factor for parameters specified by r_U
0027 % Um1q     = I * d_I^2 diagonal blocks of \bar U_1 (numerator of muv1^2)
0028 % Um2      = I * d_I^2 diagonal blocks of U2
0029 %
0030 %
0031 % Wolfgang Förstner
0032 % wfoerstn@uni-bonn.de
0033 %
0034 % See also GaussHelmertModel
0035 
0036 function [vgi,Xv,Rgm,nabla_cg,muv,muv1,U1mqc,U2mc]=diagnostics_GHM_constraints_multi_d...
0037     (r_U,I,d_I,d_G,Am,Bm,Cov_xx,W_xx,W_gg,vv)
0038 
0039 [~,U] = size(Am);
0040 
0041 delta_0=4.13;  % all sizes referring to d_J.
0042 
0043 if ~isempty(r_U)
0044     % ranges
0045     r_C= r_U;
0046     r_D= setdiff(1:U,r_U);
0047     % partitioned design matrix A=[C,D]
0048     Cm = Am(:,r_C);
0049     Dm = Am(:,r_D);
0050     % reduced design matrix
0051     Cmr    = Cm - Dm / W_xx(r_D,r_D) * W_xx(r_D,r_C);  % C_reduced
0052     Cov_11 = Cov_xx(r_C,r_C);
0053 end
0054 
0055 %% observational groups
0056 Rgm   = zeros(I,d_G^2);
0057 U1mqc = zeros(I,d_G^2);
0058 U2mc  = zeros(I,d_G^2);
0059 Xv    = zeros(I,1);
0060 vgi   = zeros(I,d_G);
0061 nabla_cg = zeros(I,1);
0062 muv   = zeros(I,1);
0063 muv1  = zeros(I,1);
0064 
0065 for i=1:I
0066     % detectability, test statistics, sensitivity wrt all parameters
0067     i_range     = d_I*i-(d_I-1):d_I*i;
0068     g_range     = d_G*i-(d_G-1):d_G*i;
0069     Amt_i       = Am(g_range,:);
0070     Bmt_i       = Bm(i_range,g_range)';
0071     %W_gg_ii     = inv(Bmt_i * Cov_ll_ii * Bmt_i'); % Weight matrix of group
0072     W_gg_ii     = W_gg(g_range,g_range); % Weight matrix of group
0073     Cov_gg_ii   = inv(W_gg_ii);
0074     Cov_vgvg_ii = Cov_gg_ii - Amt_i * Cov_xx * Amt_i';
0075     W_vgvg_ii   = inv(Cov_vgvg_ii+eps);
0076     % redundancy matrix of constraints
0077     R_gg_ii     = full(Cov_vgvg_ii * W_gg_ii);
0078     R_gg_ii_inv = inv(R_gg_ii);
0079     Rgm(i,:)    = R_gg_ii(:)';                  % vectorized Rgg
0080     vgi(i,:)    = Bmt_i*vv(i_range);       % residual of group
0081     Xv(i)       = vgi(i,:) * W_vgvg_ii * vgi(i,:)';                         %#ok<MINV>
0082     % Xhi^2-test statistic
0083     nabla_cg(i) = delta_0 * ...
0084         sqrt(real(max(eig( R_gg_ii \ inv(W_gg_ii) ))));
0085     % minimum detectable outlier
0086     muv(i)      = sqrt(real(max(eig(R_gg_ii_inv - eye(d_G) ))));
0087     % sensitivity factors
0088     
0089     % sensitivity wrt selected parameters range rU
0090     
0091     if ~isempty(r_U)
0092         % U's for constraints
0093         U1qc       = Cmr(g_range,:) * Cov_11  * Cmr(g_range,:)' * W_gg_ii;
0094         U1mqc(i,:) = U1qc(:)';
0095         U2c        = Dm(g_range,:)  / W_xx(r_D,r_D)  * ...
0096             Dm(g_range,:)' * W_gg_ii;
0097         U2mc(i,:)  = U2c(:)';
0098         % check      = U1qc + U2c + R_gg_ii - eye(d_G)
0099         % I=U_C_bar+U_D+R, check-eye(d_I)
0100         
0101         % sensitivity factors
0102         muv1(i)   = sqrt(real(max(eig(U1qc / R_gg_ii ))));
0103     end
0104 end
0105

Generated on Sat 21-Jul-2018 20:56:10 by m2html © 2005