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

 params see GaussHelmertModel

 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 % params see GaussHelmertModel
0008 %
0009 % Wolfgang Förstner
0010 % wfoerstn@uni-bonn.de
0011 %
0012 % See also GaussHelmertModel
0013 
0014 function [vgi,Xv,Rgm,nabla_cg,muv,muv1,U1mqc,U2mc]=diagnostics_GHM_constraints_multi_d...
0015     (r_U,I,d_I,d_G,Am,Bm,Cov_xx,W_xx,W_gg,vv)
0016 
0017 [~,U] = size(Am);
0018 
0019 delta_0=4.13;  % all sizes referring to d_J.
0020 
0021 if ~isempty(r_U)
0022     % ranges
0023     r_C= r_U;
0024     r_D= setdiff(1:U,r_U);
0025     % partitioned design matrix A=[C,D]
0026     Cm = Am(:,r_C);
0027     Dm = Am(:,r_D);
0028     % reduced design matrix
0029     Cmr    = Cm - Dm / W_xx(r_D,r_D) * W_xx(r_D,r_C);  % C_reduced
0030     Cov_11 = Cov_xx(r_C,r_C);
0031 end
0032 
0033 %% observational groups
0034 Rgm   = zeros(I,d_G^2);
0035 U1mqc = zeros(I,d_G^2);
0036 U2mc  = zeros(I,d_G^2);
0037 Xv    = zeros(I,1);
0038 vgi   = zeros(I,d_G);
0039 nabla_cg = zeros(I,1);
0040 muv   = zeros(I,1);
0041 muv1  = zeros(I,1);
0042 
0043 for i=1:I
0044     % detectability, test statistics, sensitivity wrt all parameters
0045     i_range     = d_I*i-(d_I-1):d_I*i;
0046     g_range     = d_G*i-(d_G-1):d_G*i;
0047     Amt_i       = Am(g_range,:);
0048     Bmt_i       = Bm(i_range,g_range)';
0049     %W_gg_ii     = inv(Bmt_i * Cov_ll_ii * Bmt_i'); % Weight matrix of group
0050     W_gg_ii     = W_gg(g_range,g_range); % Weight matrix of group
0051     Cov_gg_ii   = inv(W_gg_ii);
0052     Cov_vgvg_ii = Cov_gg_ii - Amt_i * Cov_xx * Amt_i';
0053     W_vgvg_ii   = inv(Cov_vgvg_ii+eps);
0054     % redundancy matrix of constraints
0055     R_gg_ii     = full(Cov_vgvg_ii * W_gg_ii);
0056     R_gg_ii_inv = inv(R_gg_ii);
0057     Rgm(i,:)    = R_gg_ii(:)';                  % vectorized Rgg
0058     vgi(i,:)    = Bmt_i*vv(i_range);       % residual of group
0059     Xv(i)       = vgi(i,:) * W_vgvg_ii * vgi(i,:)';                         %#ok<MINV>
0060     % Xhi^2-test statistic
0061     nabla_cg(i) = delta_0 * ...
0062         sqrt(real(max(eig( R_gg_ii \ inv(W_gg_ii) ))));
0063     % minimum detectable outlier
0064     muv(i)      = sqrt(real(max(eig(R_gg_ii_inv - eye(d_G) ))));
0065     % sensitivity factors
0066     
0067     % sensitivity wrt selected parameters range rU
0068     
0069     if ~isempty(r_U)
0070         % U's for constraints
0071         U1qc       = Cmr(g_range,:) * Cov_11  * Cmr(g_range,:)' * W_gg_ii;
0072         U1mqc(i,:) = U1qc(:)';
0073         U2c        = Dm(g_range,:)  / W_xx(r_D,r_D)  * ...
0074             Dm(g_range,:)' * W_gg_ii;
0075         U2mc(i,:)  = U2c(:)';
0076         % check      = U1qc + U2c + R_gg_ii - eye(d_G)
0077         % I=U_C_bar+U_D+R, check-eye(d_I)
0078         
0079         % sensitivity factors
0080         muv1(i)   = sqrt(real(max(eig(U1qc / R_gg_ii ))));
0081     end
0082 end
0083

Generated on Mon 08-Jan-2018 17:21:49 by m2html © 2005