rstap: An R Package for Spatial Temporal Aggregated Predictor Models
Adam Peterson Brisa Sanchez

	
	Transition from circular-ribbon to parallel-ribbon flares associated with a bifurcated magnetic flux rope
    Z. Zhong1,2, Y. Guo1,2, M. D. Ding1,2, C. Fang1,2, Q. Hao1,2
    December 30, 2023
=========================================================================================================


    The rstap package implements Bayesian spatial temporal aggregated predictor models in R using the probabilistic programming language .  A variety of distributions and link functions are supported, allowing users to fit this extension to the generalized linear model with both independent and correlated outcomes. 
	
    
	

§ INTRODUCTION

	
     The built environment refers to the human made space in which people live, work, and recreate on a day-to-day basis <cit.>. The features of the built environment are many - ranging from sidewalk availability, street density, green space, ambient light or sound levels to the physical presence of amenities like community centers or businesses that can be mapped as point locations according to their address.  In this paper we are primarily concerned with the latter and for simplicity refer to them as built environment features (BEFs). An expanding body of research is focused on quantifying the health impact of BEFs since they constrain and enable everyday choices that may contribute to the development of disease <cit.>. For example, in the United States, a preponderance of convenience stores or fast food restaurants surrounding schools offering predominantly “junk”  food may increase the odds an attending student develops obesity
   <cit.>.A key limitation in this research area is that the spatial and temporal scales at which BEFs may impact health are unknown; that is, there is uncertainty about the domain, in either space or time, at which BEFs are most relevant for health. One method in current use that has focused on addressing questions about the spatial scale is the distributed lag model (DLM) <cit.>. DLMs use counts of businesses within a discrete series of user-specified distance radii to examine how effects of BEFs decay as a function of distance between BEFs and study participants. However, DLMs are limited by the fact that they do not specifically estimate a spatial scale parameter and require a sufficient number of distance radii parameters <cit.>. Peterson et al <cit.> recently proposed spatial-temporal aggregated predictor (STAP) models to address questions about both the spatial-temporal scales and to avoid the accuracy issues associated with discretizing distances in the DLM.  The purpose of the present article is to provide a general overview of the  package (version 1.0.2), which implements STAP models. This package utilizes Hamiltonian Markov Chain Monte Carlo (HMC) estimation via <cit.>.  can be considered as being in the same family of packages as  <cit.> and  <cit.> in the sense that it utilizes  as the estimation engine to fit models of a designated form. We begin with a review of the modeling framework followed by an introduction of the software using both simulated and real data.  Throughout we use food outlets as examples of BEFs - e.g. fast food restaurants and coffee shops - and body mass index (BMI) or obesity as examples of outcomes. We end by discussing current limitations and future plans for extending the package and modeling framework.
	

§ MODEL DESCRIPTION

	Let Y_i be a vector containing n_i repeated outcome measures for subject i= 1, 2, ..., N, with corresponding mean vector μ_i = [μ_1,...,μ_n_i]. The modeling objective is to estimate a generalized linear mixed model (GLMM) regression of Y_i on typical covariates and summaries of built environment features. Define Z_i as a matrix of covariates with corresponding population parameters δ for said regression and let b_i ∼ MVN(0,Σ) be random effects for subject i, with corresponding design matrix W_i. The STAP model extends the standard GLMM, g(μ_i) = α + Z_i δ + W_i b_i, where g(·) is a link function, by estimating the latent effect of BEFs on μ_i. For simplicity, we describe the model focusing on one type of BEF, fast food restaurants (FFR), and univariate outcome μ_i ∈ℝ^1, discussing how it can be expanded later. In order to estimate the BEF effect and spatial scale, the STAP model requires pairwise distances d between each subject and each BEF, in addition to a weighting function  𝒦_s chosen so that K_s(0,θ) = 1, lim_d →∞𝒦_s(d,θ) = 0 where distance d ∈ [0,∞) and spatial scale θ∈ℝ^+. This corresponds to the substantive belief that a given BEF's maximum impact is made when a subject is as close as possible to it and vice versus. Let 𝒟_i be the set of all aforementioned pairwise distances, d, between all the BEF locations and subject i. The model is then:
    
    g(μ_i)    = α + β X_i(θ) + Z_i^Tδ1
    
        X_i(θ)    = ∑_d ∈𝒟_i𝒦_s ( d/θ^s )

    Hence, X_i represents subject i's  cumulative exposure to the particular type of BEF, accumulated according to the spatial scale θ^s. The coefficient β is the estimated difference in g(Y_i) associated with a one unit increase in X_i(θ) all else equal. Equivalently, β can be thought of as the change in g(Y_i) associated with placing one new BEF at a distance of 0 units from the subject, all else equal, since this placement results in one unit higher X_i(θ). A larger spatial scale corresponds to a BEF that continues to have an impact at larger distances and vice versus. One may consider having multiple types or classes of BEFs in the same model. For instance, consider FFRs in addition to convenience stores and gyms. The model is extended by indexing the classes of BEFs by j=1, …, J, and re-writing  model (1) as: 
     
    g(μ_i|b_i)    = α + ∑_j=1^J β_j X_i,j(θ_j)+ Z_i^Tδ  + W_i^Tb_i 2
    
        X_i,j(θ_j)    = ∑_d ∈𝒟_ij𝒦_s ( d/θ^s_j )

      where 𝒟_ij is the set of distances from individual i to BEFs of class j associated with subject i, θ_j is a spatial scale specific to the jth class and g(μ_i|b_i) is the transformed mean of outcome Y_i conditional on b_i, g(E[Y|b_i]). Extending this framework to a setting in which μ_i ∈ℝ^n_i and defining a model with both temporal and spatial components under the STAP framework yields the following expression:
    
    g(μ_i|b_i)    =  α + ∑_j=1^J β_j X_i,j(θ)  +Z_iδ + W_ib_i 3
    
        X_i,j,k(θ)    = ∑_(d,t) ∈𝒟_ijk𝒦_s ( d/θ^s_j )𝒦_t(t/θ^t_j)

    Here t ∈ [0,∞) represents the time subject i has spent “exposed” to BEFs of class j at distance d by the kth measurement visit (k = 1,..., n_i). Thus, X_i,j,k(θ) represents the exposure of the ith subject accumulated by the kth visit attributable to all BEFs in the jth class assuming the spatial and temporal components of the BEF exposure are independent.  The scales for BEF of class j, are θ_j = (θ_j^s,θ_j^t) where θ_j^t represents the temporal scale at which exposure to a given BEF is associated with g(μ_i).  In contrast to 𝒦_s, the weight function 𝒦_t is chosen so that 𝒦_t(0)=0,lim_t→∞𝒦_t(t)=1; β_j in this setting is interpreted as the difference in g(μ_i) when a BEF is placed at distance 0 from the subject, for an amount of time that approaches infinity. Note here that the statement regarding infinite time will always be true but the BEF may exert 99% of its impact in a finite amount of time when the temporal scale θ_j^t is small. 
    
    

§ PARAMETER ESTIMATION

    
    All the above models are fit according to a Bayesian framework in the  package, allowing the incorporation of prior information into the model and the final goal being posterior inference. To that end, the  package utilizes  to fit models, drawing samples from the posterior distribution of model parameters obtained from an extension of the static HMC sampler <cit.>, the Generalized No-U Turn Sampler (NUTS) <cit.>. “Hamiltonian Monte Carlo (HMC) is a Markov chain Monte Carlo (MCMC) algorithm that
avoids the random walk behavior and sensitivity to correlated parameters that plague many
MCMC methods by taking a series of steps informed by first-order gradient information.
These features allow it to converge to high-dimensional target distributions much more
quickly than simpler methods such as random walk Metropolis<cit.> or Gibbs sampling” <cit.>. This trade off comes at the expense of calculating the gradient of parameters to be estimated and subsequent numerical integration to solve the Hamiltonian differential equations <cit.>. Consequently this framework both enables and constrains the kinds of models that can currently be fit in . Since each X_i,j is a function of θ_j, several distances and/or times, the time required to differentiate θ_j and β_j via 's automatic reverse mode differentiation scheme is longer than it would be in a typical regression setting. The computational complexity of this differentiation is a function of both the specification of 𝒦_s and/or 𝒦_t, as well as the number of BEFs.  To provide a bound on this complexity in , the number of BEFs is constrained to those that lie within a pre-specified maximum distance. This maximum distance also enforces an upper bound on the possible values of θ^s. Since the monotonic function 𝒦_s discussed previously would become uniform for an infinitely large scale θ^s evaluated at finite distance, an upper bound must be chosen to ensure the posterior of the scales do not tend toward infinity in the case where the estimated spatial exposure approaches uniform on the domain of distances under which the model was fitted. The current setting in  divides the max distance by the maximum ninety seventh quantile of the chosen 𝒦 function to obtain this upper bound - an explicit formulation of this bound can be found in the Appendix, Figure <ref>. This allows users to recognize the “approximate” uniformity of the exposure effect of a BEF and re-adjust the inclusion distance as needed. An example of this is given on the  website (https://biostatistics4socialimpact.github.io/rstap/articles/Introduction.html).With our simulated and real datasets examining between one and five STAPs with anywhere from a few hundred to a few thousand subjects, and anywhere from a few dozen to a few hundred pairwise distances per person,  sampling times for 500 posterior samples - including warmup - vary from a few seconds to just under an hour. Further comments on these constraints and future work for plans to speed these sampling times are elaborated upon in the discussion. Similar to  and ,  also offers draws from the posterior predictive distribution and pointwise log-likelihood. The former allows assessment of model fit between the model's outcome prediction's, Ŷ, and Y, while the latter permits model selection via Watanabe-Akaike information criteria <cit.>. Software functions in  that take advantage of these features are discussed in our example analyses.
    

§ SOFTWARE

     The current development version of  can be downloaded and installed from GitHub via  . Models are fit in  using the following procedure, which is also summarized in Figure 1. This framework is akin to the  and  functions in the popular  package <cit.> as well as the  function in  for correlated data <cit.>. A user provides a formula which specifies the kind of model, the family of distributions and the link function under which the mean of the distribution is hypothesized to be related to the linear predictors. Currently Binomial, Bernoulli, Poisson and Gaussian distribution families are supported by .In contrast to  and , a total of two or three data sets must be supplied to the function call, each containing an ID to relate them to each other.  One of these datasets, “subject_data”, must always be supplied to the function with a unique subject identifier on each row, any potential group IDs, and any other standard covariates to include in the model.  The other two datasets  contain the built environment information: specifically, the distances or times between subjects and BEFs. One or both of these are required dependent upon whether only or both spatial and temporal components are specified in the model. The typical structure for built environment datasets - “distance_data” and “time_data” - is given below in Table <ref>.  Each row defines a unique subject BEF association with its corresponding distance. The  argument in this model would be the, identically named, “subject_ID” string, since that is the name of the column containing the unique subject ID that maps subjects' distances in the “distance_data” data frame to their appropriate  covariate information in the “subject_data” data frame. In the case of grouped data, an additional column would define the group to which the data should be associated - see Tables <ref> and <ref> for examples.These data are then passed to one of several pre-set compiled  programs akin to  for sampling. After this a  object will be defined using the samples. The  class is modeled off of the  class in  allowing a user to print out a clean summary of the model estimates including diagnostics of the MCMC fit, such as MCMC standard errors and standard split-chain R̂ diagnostics <cit.>.

    

 §.§ A Worked Example

    There are two main functions in the  package:  and . Two other functions included in the package,  and  have the same purpose as the previous, but simply set the family option to be , implying a gaussian distribution with an identity link will be used to fit the data. Suppose one has a subject level dataset with subjects' sex and BMI recorded, as well as a dataset with the distances, in miles, between each subject and all FFRs within an area of, say, 5 miles of each subject. Then the code in Figure <ref> would correspond to the following model and prior distribution specifications, relating the average BMI to subjects' sex and exposure to FFRs:
    
    BMI_i    = α + δSex_i + βFast_Food_i(θ^s) + ϵ_i
    ϵ_i    iid∼ N(0,σ)     σ∼ C^+(0,5)
    β   ∼ N(0,4)     δ∼ N(0,4)
    α   ∼ N(26,4)     log(θ^s) ∼  N(1,1)

    Further details regarding the specification of the priors are discussed in Section 4.3.
    
    
    

    

 §.§ Formula

    
     The formula argument to all  functions contains the relationship between the response and the predictors, as well as the specific configuration for the spatial-temporal predictors or group terms included in the model. Group terms are specified identically to 's syntax, where  denotes that a random intercept and slope for “coefs” at the level of “group” should be included in the model. Spatial-Temporal components may be specified in a number of ways: Table <ref> below shows the different syntax that may be used to fit a model akin to the one in Section 4.1, with differing spatial-temporal components using the two different weight functions implemented in . Note that while ,  and  each refer to different kinds of STAPs, references made to s hereafter encompass all three kinds of covariates unless stated otherwise.  



  
    

 §.§ Prior distributions

    
    Discussion of how priors should be placed on parameters in standard regression models has already been elaborated in the packages mentioned previously <cit.>. In the next two sections we briefly discuss new considerations and  syntax for setting priors on the STAP specific variables: θ^s, θ^t, β.
    
    

  §.§.§ Spatial-Temporal Scale Parameters

     allows for custom specification of the priors on any scale parameters, θ^s,θ^t . This allows each scale to have its own prior distribution according to the user's desire. The syntax used in Figure <ref>, corresponds to assigning the same prior distribution to all spatial-temporal scales, θ^s,θ^t. In cases where there is only one spatial predictor or there are similar a priori beliefs about the spatial-temporal scales, this allows for easy assignment of scale priors. Alternatively, if one were interested in assigning different priors, perhaps when fitting a model with both spatial and temporal components on the class of “Fast Food” restaurants from the previous example, the model and corresponding code would look as follows:
    
    BMI_i    = α + δSex_i  + βFast_Food_i(θ^s,θ^t)  + ϵ_i
    log(θ^s)    ∼ N(1,1)     log(θ^t) ∼ N(1,2)

    
    Currently only log and folded normal prior distributions are supported for θ^s,θ^t in  as these are the only priors with which we have done extensive simulations and testing. These can be used for less and more informative priors, respectively.
    


  §.§.§ Population and group level parameters

 Population and group level parameters are constructed and assigned priors in essentially the same form as they are in the  <cit.> package. The only difference is that there are separate  and  arguments to allow for differing prior specifications for the δ and β parameters, respectively. Note that while the user can choose to autoscale priors for the standard regression coefficients as in , priors are always set on standardized scales for β parameters. That is, during estimation, the estimated X(θ) values are mean centered and scaled to allow for easier computation, making the parameter space for a β's far easier to traverse for the sampler <cit.>. Further, this allows placement of equivalent priors to be set on multiple  covariates that might otherwise be inappropriate due to highly differing exposure levels <cit.>. 

 
    
    

§ CONTROL

    Regular control of the sampler will occur via three arguments: , , and . The first controls the total number of iterations for which the sampler is run. The second controls the subset of those previous iterations that are used to “warm-up” or tune the number of steps and mass-matrix which are used to propose new samples via NUTS<cit.>. Finally,  controls the target sample acceptance ratio which determines the resolution at which the sampler explores the posterior. A higher  corresponds to a higher target ratio and consequently, a longer sampling time. Two additional arguments,  and , specifies the number of processors to use and the number of chains to run. For example, if  and  then two chains will be run in parallel via the  package <cit.>. Further specification of 's algorithm specific variables occurs via the optional argument  analogous to  <cit.>.
    
    

§ EXAMPLE ANALYSES

    In the following four sections, we demonstrate the STAP framework in settings involving continuous and binomial outcomes from both simulated and real data.

    

 §.§ Simulation I: Uncorrelated Outcomes

    We demonstrate the  package first on a simulated dataset where the built environment features are homogeneously distributed in a 3 x 3 square, and the subjects are similarly homogeneously distributed in a 1 x 1 square centrally interior to the latter. Although obviously fictitious, this set-up provides some intuition as to how the different components of the model can be set-up and interpreted. The spatial arrangement and distribution of pairwise distances between subjects and BEFs can be seen in Figure <ref>. 

    
    The outcome is then simulated according to the following model, with the intention of mimicking the marginal distribution of BMI in a random sample of subjects. Labeling our simulated BEFs as FFRs and simulating the spatial exposure using the complementary error function,  we recreate the example previously discussed in 4.1. Visualizations of the spatial decay and resulting exposure distribution can be seen in Figure <ref>.
    
    BMI_i    = 22.5 - Sex_i 0.8 + Fast_Food(θ^s=.5) 1.2 + ϵ_i
    ϵ_i    ∼ N(0,2.3)

    

    
    
    Assuming the data has been appropriately structured as described in Section 4, the model can then be fit with the syntax in Figure <ref> after loading the  library. Typical model printout is seen below in Figure <ref>.
    
    Calling  on the fit object will produce a longer print out, with relevant convergence diagnostics and WAIC, if specified - Appendix Figure <ref>. Furthermore, typical model fit diagnostics such as posterior predictive checks <cit.> can be accessed directly through the   function and the  function from the  package<cit.>. The posterior predictive checks and NUTS energy diagnostics <cit.> - available again through the  package - for this model can be seen in Figure 5 and the corresponding code can be found in the Appendix - Figure <ref> 

    
    

    





 §.§ Simulation II: Correlated Outcomes


Our second simulation showcases the typical STAP data structure in a longitudinal setting. Consider a model similar to the previous except we now have multiple measurements on each subject, in addition to the time each subject spent at their respective locations.

    BMI_i,j|b_i   = α + Z_i,j^Tδ + X_i,j(θ^s=.8,θ^t=18)β + b_i + ϵ_ij   j = 1,2
    
     b_i    ∼ N(0,1.5)
    ϵ_ij   ∼ N(0,2)

Prior to model fitting, the data collected for subjects and BEFs can be quite complex, as subjects move across space and time and businesses will open and close rendering a dynamic covariate space. Structurally, however, the data submitted to  will be similar to the uncorrelated spatial setting demonstrated in Simulation I, with the addition of at least one new column for a “group ID” so that distances or times for specific BEFs are associated with the subject at the appropriate measurement or group level. Note that multiple groups may be included in the time or distance data frame as coded ID column variables and, consequently, multiple IDs passed to the  “group_ID” argument. See Tables <ref> and <ref> for how the two data frames supplied to  and , respectively, are structured in a possible setting where a subject moves locations between two consecutive measurements. 
 







Simulating data under the model above with the same spatial configuration as the first and setting β = 1, we fit the model drawing two thousand samples on four independent chains. The first thousand samples from each chain are used for warm-up resulting in four thousand total posterior samples.  We use similar priors as before and place priors on the subject specific variance parameter in accordance with 's recommendations for hierarchical covariance matrices <cit.>.  The model output printed from  can be found below along with a visualization of both the spatial temporal exposure estimates in Figure <ref>. A full workup for how this kind of data may be generated and compiled can be found at the package's website along with other, more complicated, simulated spatial structures and weight functions  (https://biostatistics4socialimpact/rstap).






    

 §.§ California FitnessGram and NETS data

    We use FitnessGram data for 5th grade students in California during the 2010-2011 Academic year to examine associations between availability of fast food restaurants near their school and children's obesity. Publicly available data files from the California Department of Education (CDE) contain information on students' obesity status, grouped at the school level, or within sex and race/ethnicity groups within the school<cit.>. We also obtained school-level covariates from the CDE, namely, charter status, percent of children eligible for free or reduce priced meals, and majority race/ethnicity of the school. We joined this data set to the Census track characteristics where the school was located, to include school-neigbhorhood characteristics: median household income and level of urbanization. Level of urbanicity around the school was coded as Urban, Suburban and Rural. For this analysis, we subset the data to only examine schools within urban areas.
    
    Participating in the FitnessGram test is required by the State of California, and as such, informed consent is not required. Moreover, since all personal identifiers are removed by the CDE prior to making the school-level aggregated data publicly available to researchers, this research is considered exempt from ethics review by the institutional review boards of The University of Michigan.
    
    The locations of California FFRs that were open during this time period were obtained from a commercial source<cit.>. Distances were calculated between each school and businesses within ten miles categorized as “Fast-Food Chain” or “Fast-Food Non-Chain”. These two business classes were combined to form one class, “Fast Food Restaurants“ (FFR), to be included in the model. The ten mile inclusion distance was chosen as a conservative upper bound as previous associations estimated with a DLM decline to zero at far shorter distances <cit.>.
    
    We fit the following model in order to estimate the odds of obesity for a fifth grade child  as function of their Fast Food environment, after adjusting for the aforementioned confounders:
    
    
    logit(μ_i)    = α + Z_i^Tδ + X_i,FFR(θ_s)β
    α   ∼ N(0,3)      δ∼ N(0,2) 
    β   ∼ N(0,2)     log(θ^s) ∼ N(1,1)

     In this example we chose priors for the prior intercept and regression coefficients according to Gelman et al's principles for weakly informative priors in a logistic regression setting<cit.>. We set the prior for the STAP regression coefficient similarly using previous research <cit.> for further justification. Corresponding to our inclusion distance of 10 miles, a prior distribution was set on the spatial scale such that the median exposure was 0.01 - effectively negligible - after 5 miles but whose tail was long enough that a true higher scale would not be as severely penalized as it would with, say, a normal prior. Note that we do not use uninformative priors in this setting as it not only makes sampling times prohibitive, but corresponds to placing substantial probability mass on many scales that would set equivalent exposures to BEFs both close and extremely far away from schools.  We code this prior information and model in , with the function call seen in Figure <ref>.
    
    A visualization for the differing model's estimates of the spatial scale can be found in Figure <ref> along with the initial prior distribution and FFR effect. 


    
    

    
    We found that higher exposure to FFRs near schools is associated with a higher odds of obesity after adjusting for confounders.  In urban areas, one FFR placed zero miles from a school is associated with an odds ratio of 1.001558 (95 % CI 1.000471 to 1.002783). The small magnitude of this estimate reflects the saturation of FFRs in urban environments, i.e. one more FFR has negligible association with obesity. To emphasize this fact, the standard deviation of the X(θ^s) covariate in urban areas is 23.8 (12.59 to 51.8), implying that if 24 stores were added at zero miles from a school - one standardized effect - the odds ratio for obesity would be roughly 4% as seen in Figure <ref> (a). 

    Moreover, using our spatial scale estimates we can see that the effect of one additional restaurant is effectively negligible after, roughly,  six miles. We calculate this precisely for Table <ref>, using the  function in . This function requires a specification of what exposure to consider “effectively negligible". For Table <ref> we use a value of 0.01  (code in Appendix, Figure <ref>). That is, we find the value d for which 𝒦(d/θ̂^s) = 0.01 for the corresponding θ̂^̂ŝ estimates. 
 



 §.§ FitnessGram data - Correlated outcomes

To demonstrate how STAPs can be used with correlated binomial outcomes, we use the number of obese boys and girls at each school in the same setting as before. In terms of data structure, since the regression design matrix must be broken out into two groups, the distance data submitted to  must also be broken out, or in this case, copied into two groups to reflect the corresponding spatial temporal exposure for the observation at the appropriate group level. The data are copied in this setting because the BEF exposure is the same for the school, regardless of whether the boy or girl outcome is of interest. In this case, the modification is reflected in the  ID key and covariate shown in the  function below. This function fits the model with a random intercept and similar priors placed on τ as in our longitudinal simulation and δ as in the previous example. Additionally, the same visualizations are provided for this new model in Figure <ref>.
    
    
    logit(μ_i,j|b_i )    = α + Z_i,j^Tδ + X_i,FFR^T(θ_s)β + b_i  
    
        j     denotes an outcome for girls or boys at the  i th school
    
        b_i    ∼ N(0,τ) priors in function call

    
    
    

    
    

    
    Examining this model's estimates we can see that a similar median exposure scale is attained, with a much wider uncertainty. This is likely due to the increased number of parameters being estimated in this new model formulation, in addition to the increased uncertainty imposed by the introduction of an additional latent parameter into the model.
    
    

§ DISCUSSION

    
	 The present paper is meant to provide a general overview of the R package  implementing STAPs using the probabilistic programming language  for full Bayesian inference. This article provides a starting point for the reader interested in understanding how to fit basic STAP models.The current extension planned for the next release of  is a custom HMC sampler with symbolic derivatives of the log probability hard coded for typical STAP model configurations that, not requiring automatic differentiation, will sample from the posterior much more quickly. This will be advantageous for modeling scenarios in which there are a large number of STAPS to be fit and/or a large number of built environment features associated with each subject. Other extensions could include providing functions for modeling stap-interactions in a setting where it is believed a scale may differ among levels of some categorical, e.g. gender, covariate.
	  

§ ACKNOWLEDGEMENTS

	  The authors would like to acknowledge funding sources RO1 HL131610 and RO1 HL136718 in support of this work.

 

unsrt


 

§ APPENDIX




 §.§ Spatial Temporal Scale Upper Bound Construction

    For weighting function 𝒦_j(d) associated with each STAP j=1,...,J we find θ^up, the upper bound used in estimating θ_j^s, j = 1,...,J. Here, d_j^* is the distance specific to the specified 𝒦 for class j at which the evaluation is equal to 0.975. Since different staps may have different weight functions, and consequently, different d_j^*'s taking the minimum across these distances results in a unified upper bound for all θ, which is convenient computationally.
    









