[ncl-talk] Spatial plot for Significance of Linear Trend

Sri Nandini snandini at marum.de
Mon Mar 5 02:36:57 MST 2018


Hello

I am trying to calculate statistical significance of a linear trend (Regression analysis).
The regression analysis simple enough.

I have troubles understanding how the probabilities are computed based on the trend i calculate?
What i actually wanted to do is calculate the 2 tailed ttest probabilities but it doesn't seem possible with NCL, unless someone else knows of a way?
Or does using the NCL function "betainc" to get the p-value
for a student test is the same thing?


>>>My code is below:


f = addfile ("air.mon.mean.nc", "r")
  latS   =  29.3
  latN   =  67.26 
  lonL   = 19.
  lonR   =  60.

  TIME1   = f->time
  YYYY1   = cd_calendar(TIME1,-1)/100                 ; entire file
  iYYYY1  = ind(YYYY1.ge.yrStrt .and. YYYY1.le.yrLast)
  t2m    = f->air(iYYYY1,:,:)
  printVarSummary(t2m)

  T2M    = month_to_season (t2m, season)
  printVarSummary(T2M)
  nyrs   = dimsizes(T2M&time)

  t2mw     = T2M({lat|latS:latN},{lon|lonL:lonR},time|:)
  printVarSummary(t2mw)

  t2mw = dtrend(t2mw,False)
  printVarSummary(t2mw)

  eofT = eofunc_Wrap(t2mw, neof, optEOF)
  printVarSummary(eofT)

;==========================================================================
eof_regres= eofT                            ;create an array with meta data 
do ne=0,neof-1
eof_regres(ne,:,:)=  regCoef(eof_ts(ne,:), t2mw(lat|:,lon|:,time|:));
  end do
  printVarSummary(eof_regres)

;==========================================================================
;Calculate significance at 95% confidence
;The regCoef function returns the following attributes: yintercept (y-intercept), tval (t-statistic), rstd (standard error of the estimated regression coefficient) and nptxy (number of elements used). 
;==========================================================================

    df = onedtond(eof_regres at nptxy,dimsizes(eof_regres))-2        ; degrees of freedom
    tval = onedtond(eof_regres at tval, dimsizes(eof_regres))        ; t-statistic

    b = tval                                  ; b must be same size as tval (and df)
    b = 0.5
    prob = betainc(df/(df+tval^2),df/2.0,b)
    printVarSummary(prob)                     ; prob(nlat,nlon)

   ;copy_VarCoords(eof_regres,prob)
   prob at long_name = "probability"
   ;prob at units = "fraction: [0,1]"

   prob_95=where (prob.ge.0.95,1,0)           ; **** note 0.95 ****
   copy_VarCoords(prob, prob_95)
   printVarSummary(prob_95) 

   eof_regres at long_name   = "regression coefficient"
   prob at long_name = "probability"
  
;==================================The error message i get is:

warning:VarVarWrite: rhs has no dimension name or coordinate variable, deleting name of lhs dimension number (1) and destroying coordinate var,  use "(/../)" if this is not desired outcome
warning:VarVarWrite: rhs has no dimension name or coordinate variable, deleting name of lhs dimension number (2) and destroying coordinate var,  use "(/../)" if this is not desired outcome
warning:["Execute.c":8640]:Execute: Error occurred at or near line 135 in file NAOregreT2M_obs.ncl


Variable: eof_regres
Type: float
Total Size: 49856 bytes
            12464 values
Number of Dimensions: 3
Dimensions and sizes:    [evn | 2] x [76] x [82]
Coordinates: 
            evn: [1..2]
Number Of Attributes: 11
  tval :    <ARRAY of 6232 elements>
  yintercept :    <ARRAY of 6232 elements>
  rstd :    <ARRAY of 6232 elements>
  nptxy :    <ARRAY of 6232 elements>
  long_name :    EOF: DJF: Monthly mean of surface temperature
  _FillValue :    -9.96921e+36
  method :    transpose
  matrix :    covariance
  pcvar :    ( 59.16475, 20.32462 )
  eval :    ( 18982.65, 6521.03 )
  eval_transpose :    ( 169.488, 58.22348 )



-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180305/2d0554db/attachment.html>


More information about the ncl-talk mailing list