Announcement

Collapse
No announcement yet.
X
  • Filter
  • Time
  • Show
Clear All
new posts

  • Numerical derivatives with deriv() when using factor variables

    Dear all,
    We are writing a program that involves calculating second derivatives numerically using mata's deriv() function. We want to allow factor variables to enter the function for which the second derivatives are calculated. Creating a mata matrix from stata factor variables means that this matrix will contain columns of zeros. In this case, deriv() aborts with the error message "could not calculate numerical derivatives -- flat or discontinuous region encountered", which makes perfect sense, since the function needs to be flat with respect to the coefficients of the omitted (base level) variables. See the code below for an illustration of this problem using a super-simple probit example..

    Code:
    clear
    set seed 152507
    set obs 99
    gen byte x = rbinomial(1,0.5)
    gen byte y = rbinomial(1,normal(-0.5+x))
    gen byte one = 1
    
    mata : XF = st_data(.,"i.x one")
    mata : X = st_data(.,"x one")
    mata : Y = st_data(.,"y")
    
    mata
    : void llprobit_v(c, RHS, LHS, lnf)
    {
     lnf = ln(normal(RHS*c')) :* (LHS) :+ ln(normal(-RHS*c')) :* (1 :- LHS)
     }
    end
    
    probit y i.x
    
    mata
     D = deriv_init()
     deriv_init_evaluator(D, &llprobit_v())
     deriv_init_evaluatortype(D, "v")
     deriv_init_params(D, st_matrix("e(b)"))
     deriv_init_argument(D, 1, XF)  
     deriv_init_argument(D, 2, Y)  
     deriv(D, 2)
    end
    However, one can make deriv() do what one want it to do, i.e. set the (second) derivatives with respect to the coefficients of the omitted variables to zero, by specifying deriv_init_search(D, "off"), i.e. preventing deriv() from searching for optimal delta values. In the "playing around" example above, this works well, i.e. the derivatives of the other coefficients are the same as those obtained by not specifying x as a factor variable and not specifying deriv_init_search(D, "off").

    Code:
    mata
     D = deriv_init()
     deriv_init_evaluator(D, &llprobit_v())
     deriv_init_evaluatortype(D, "v")
     deriv_init_search(D, "off")
     deriv_init_params(D, st_matrix("e(b)"))
     deriv_init_argument(D, 1, XF)  
     deriv_init_argument(D, 2, Y)  
     deriv(D, 2)
    end
    However – and not surprisingly - in more complex applications, specifying deriv_init_search(D, "off") seems to cause deriv() to calculate some rather 'strange' values for the derivatives.
    So my question is:is there another way to get deriv()to set the derivatives with respect to the coefficients of omitted variable coefficients to zero?
    (This question is obviously related to the problem of how to deal with factor variables when using optimize(). However, the solution of 'telling' optimize() by specifying optimize_init_constraints()that the respective coefficients are constrained to zero cannot, as far as I got it, be applied to deriv()).

    Thank you all for your suggestions,
    Harald
    Last edited by Harald Tauchmann; 15 Nov 2024, 09:14.

  • #2
    Thank you for providing a working example.

    deriv() does not handle constraints. optimize() and moptimize() provide their own evaluators to deriv(), and it is these evaluators that handle the constraints for deriv().

    In the following, I show the changes one might employ to get deriv() to compute the second order numerical derivatives. The changes involve applying the contraints to the predictors matrix and point estimates that are supplied to deriv(), then reversing the transformations on the resulting derivatives matrix.

    My changes take advantage of the fact that probit posts the constraints implied by the factor variables to e(Cns). Command matcprob uses this matrix to construct the tranformation matrix T and offset vector a. For more information on this, in Stata type help makecns.

    Code:
    clear all
    
    set seed 152507
    set obs 99
    gen byte x = rbinomial(1,0.5)
    gen byte y = rbinomial(1,normal(-0.5+x))
    gen byte one = 1
    
    probit y i.x
    
    * get the matrices that apply the linear constraints
    matcproc T a C
    matrix list T
    matrix list a
    matrix list C
    
    mata:
    
    // evaluator is not changed
    void llprobit_v(c, RHS, LHS, lnf)
    {
        lnf = ln(normal(RHS*c')) :* (LHS) :+ ln(normal(-RHS*c')) :* (1 :- LHS)
    }
    
    XF = st_data(.,"i.x one")
    X = st_data(.,"x one")
    Y = st_data(.,"y")
    
    // Mata copy of linear constraint tranformation matrices
    T = st_matrix("T")
    a = st_matrix("a")
    
    // values at which to compute the numerical derivatives
    b0 = st_matrix("e(b)") 
    
    // transform the data and values at which to compute the numerical derivatives
    tXF = XF*T
    tb0 = (b0 :- a)*T
    
    D = deriv_init()
    deriv_init_evaluator(D, &llprobit_v())
    deriv_init_evaluatortype(D, "v")
    // supply transformed data
    deriv_init_argument(D, 1, tXF)  
    deriv_init_argument(D, 2, Y)  
    
    // supply transformed values
    deriv_init_params(D, tb0)
    
    // compute second order derivatives
    td2 = deriv(D, 2)
    // show
    td2
    // invert
    itd2 = invsym(-td2)
    
    // transform back to unconstrained space
    V = T*itd2*T'
    V
    
    // show that above are equivalent to the VCE from -probit-
    st_matrix("e(V)")
    
    end

    Comment


    • #3
      Dear Jeff,
      thank you so much. This helps us a lot.
      Best wishes,
      Harald

      Comment

      Working...
      X