TABLE OF CONTENTS

CinitRoutine/addInitSmoothnessPenaltyGr [ Functions ]

NAME

    addInitSmoothnessPenaltyGr --- adds gradient of smoothness penalty to an input vector

FUNCTION

Computes the gradient of the smoothness penaly over the spline parameters and adds it to an input vector.

SYNOPSIS

540 void addInitSmoothnessPenaltyGr(double *gr, double * penpar, double *P, int *penaltyType,
541         double * sigma2, int *nj)

INPUTS

    gr        gradient vector
    penpar    parameters of the spline
    P         penalty matrix
    penaltyType   integer indexing the type of penalty (0=gaussian, 1=2diff, 2=2deriv, 3=log2deriv)
    sigma2    prior variance
    nj        number of B-spline basis functions

OUTPUTS

SOURCE

545 {
546     int c1 = 1;
547     if(*penaltyType==0){
548         double scalefactor = -1 / *sigma2;
549         F77_CALL(daxpy)(nj, &scalefactor, penpar, &c1, gr, &c1);
550     }
551     // second differences or second derivative
552     if(*penaltyType==1 || *penaltyType==2 || *penaltyType==3){   
553         double * temp = (double *) calloc((*nj) , sizeof(double));
554         double c1d = 1.0;
555         char uplo = 'u';
556         // compute smoothness penaly gradient.
557         F77_CALL(dsymv)( &uplo, nj, &c1d, P, nj, penpar, &c1, &c1d, temp, &c1);
558         for(int i=0; i<*nj; i++) temp[i]=temp[i]*penpar[i] / *sigma2;
559         double scalefactor=-1;
560         if(*penaltyType==3) {
561             double pen=0;
562             int c2 = 2;
563             InitSmoothnessPenalty( &pen, penpar, P, &c2, sigma2, nj);
564             scalefactor= -1/((2* *sigma2) * pen + 1);
565         }
566         // add to gradient
567         F77_CALL(daxpy)(nj, &scalefactor, temp, &c1, gr, &c1);
568         free( temp );
569     }
570 }