--- a +++ b/partyMod/src/Convenience.c @@ -0,0 +1,193 @@ + +/** + Some convenience functions + *\file Convenience.c + *\author $Author$ + *\date $Date$ +*/ + +#include "party.h" + + +/** + Linear statistic of x, y, and weights + and its conditional expectation and covariance \n + *\param x values of the transformation + *\param p dimension of the transformation + *\param y values of the influence function + *\param q dimension of the influence function + *\param weights case weights + *\param n number of observations + *\param cexpcovinf logical: recompute exp and cov of the influence fct + *\param expcovinf an object of class `ExpectCovarInfluence' + *\param ans return value; an object of class `LinStatExpectCovar' +*/ + +void C_LinStatExpCov(const double *x, const int p, + const double *y, const int q, + const double *weights, const int n, + const int cexpcovinf, SEXP expcovinf, SEXP ans) { + + C_LinearStatistic(x, p, y, q, weights, n, + REAL(GET_SLOT(ans, PL2_linearstatisticSym))); + if (cexpcovinf) + C_ExpectCovarInfluence(y, q, weights, n, expcovinf); + C_ExpectCovarLinearStatistic(x, p, y, q, weights, n, + expcovinf, ans); +} + + +/** + Moore-Penrose inverse of the covariance matrix \n + *\param linexpcov an object of class `LinStatExpectCovarMPinv' + *\param tol tolerance +*/ + +void C_LinStatExpCovMPinv(SEXP linexpcov, double tol) { + + int pq, pqn; + + /* correct working dimension */ + pq = get_dimension(linexpcov); + + /* reduce dimension to non-zero variance part */ + C_linexpcovReduce(linexpcov); + + /* reduced dimension */ + pqn = get_dimension(linexpcov); + INTEGER(GET_SLOT(GET_SLOT(linexpcov, PL2_svdmemSym), PL2_pSym))[0] = pqn; + + /* compute MPinv in reduced dimension */ + C_MPinv(GET_SLOT(linexpcov, PL2_covarianceSym), tol, + GET_SLOT(linexpcov, PL2_svdmemSym), linexpcov); + + /* make sure to reset svdmem to original dimension; + the dimension of linexpcov is reset in C_TestStatistic */ + INTEGER(GET_SLOT(GET_SLOT(linexpcov, PL2_svdmemSym), PL2_pSym))[0] = pq; +} + + +/** + Compute test statistic + *\param linexpcov an object of class `LinStatExpectCovar' + *\param type integer, 1 (maxabs) or 2 (quadform) + *\param tol tolerance +*/ + +double C_TestStatistic(const SEXP linexpcov, const int type, const double tol) { + + int pq; + double ans = 0.0; + + /* this is possibly the reduced one, see C_LinStatExpCovMPinv */ + pq = get_dimension(linexpcov); + + switch(type) { + /* maxabs-type test statistic */ + case 1: + ans = C_maxabsTestStatistic( + REAL(GET_SLOT(linexpcov, PL2_linearstatisticSym)), + REAL(GET_SLOT(linexpcov, PL2_expectationSym)), + REAL(GET_SLOT(linexpcov, PL2_covarianceSym)), + pq, tol); + break; + /* quadform-type test statistic */ + case 2: + ans = C_quadformTestStatistic( + REAL(GET_SLOT(linexpcov, PL2_linearstatisticSym)), + REAL(GET_SLOT(linexpcov, PL2_expectationSym)), + REAL(GET_SLOT(linexpcov, PL2_MPinvSym)), pq); + break; + default: error("C_TestStatistic: undefined value for type argument"); + } + + /* dimension was potentially reduced; reset to initial value */ + INTEGER(GET_SLOT(linexpcov, PL2_dimensionSym))[0] = + LENGTH(GET_SLOT(linexpcov, PL2_linearstatisticSym)); + + return(ans); +} + + +/** + Compute asymptotic conditional P-value + *\param tstat test statistic + *\param linexpcov an object of class `LinStatExpectCovar' + *\param type integer, 1 (maxabs) or 2 (quadform) + *\param tol tolerance + *\param maxpts argument to C_maxabsConditionalPvalue + *\param releps argument to C_maxabsConditionalPvalue + *\param abseps argument to C_maxabsConditionalPvalue +*/ + +double C_ConditionalPvalue(const double tstat, SEXP linexpcov, + const int type, double tol, + int *maxpts, double *releps, double *abseps) { + + int pq; + double ans = 1.0; + + pq = get_dimension(linexpcov); + + switch(type) { + /* maxabs-type test statistic */ + case MAXABS: + ans = C_maxabsConditionalPvalue(tstat, + REAL(GET_SLOT(linexpcov, PL2_covarianceSym)), + pq, maxpts, releps, abseps, &tol); + break; + /* quadform-type test statistic */ + case QUADFORM: + /* var = 0 => rank = 0 */ + if (REAL(GET_SLOT(linexpcov, PL2_rankSym))[0] > 0.5) + ans = C_quadformConditionalPvalue(tstat, + REAL(GET_SLOT(linexpcov, PL2_rankSym))[0]); + break; + default: error("C_ConditionalPvalue: undefined value for type argument"); + } + return(ans); +} + + +/** + extract the (first) response variable from a learning sample + *\param learnsample an object of class `LearningSample' +*/ + +SEXP R_get_response(SEXP learnsample) { + return(VECTOR_ELT(GET_SLOT(GET_SLOT(learnsample, PL2_responsesSym), + PL2_variablesSym), 0)); +} + + +/** + change the values of the response variable in a learning sample + *\param learnsample an object of class `LearningSample' + *\param y a REAL with new values +*/ + +void R_set_response(SEXP learnsample, SEXP y) { + + double *v, *t, *j, *dy, *p; + int i, n; + + n = LENGTH(y); + dy = REAL(y); + + if (LENGTH(R_get_response(learnsample)) != n) + error("lengths of arguments don't match"); + + v = REAL(VECTOR_ELT(GET_SLOT(GET_SLOT(learnsample, PL2_responsesSym), + PL2_variablesSym), 0)); + t = REAL(VECTOR_ELT(GET_SLOT(GET_SLOT(learnsample, PL2_responsesSym), + PL2_transformationsSym), 0)); + j = REAL(get_test_trafo(GET_SLOT(learnsample, PL2_responsesSym))); + p = REAL(get_predict_trafo(GET_SLOT(learnsample, PL2_responsesSym))); + + for (i = 0; i < n; i++) { + v[i] = dy[i]; + t[i] = dy[i]; + j[i] = dy[i]; + p[i] = dy[i]; + } +}