We performed pairwise correlation between the variables (see Supplementary Figure S1 in Supplementary Materials). Each variable consists of nt data points, where n is the number of states and t is the number of days in the target period. The two target behavioral variables, mask wearing and difference in daily social distance, were moderately negatively correlated (r = −0.41). Larger increases in the proportion of mask wearing were associated with larger decreases in distance travelled compared to the pre-COVID baseline. Both target variables were most highly correlated with the political variables (r = ±0.41–±0.53), which were themselves very highly correlated (r = 0.77–0.89). The political variables however were not very correlated with any of the mass media variables, for example, the proportion difference in people watching Fox News (r = −0.02–0.07). CDC case rate was moderately correlated with COVID- and influenza-like illnesses (r = 0.35–0.40), but not much correlated with R_t (r = −0.06). The latter might be an indication of the time lag between R_t and the onset of observable symptoms. We will re-align the data to account for this effect as identified by Granger analysis.

As mentioned, some variables might not have an immediate effect on others. We performed Granger analysis on pairs of time series variables. For each pair of variables x and y, the predictor variable x was shifted with successive time lags of 1–30 days, and the shifted previous values of x, together with the unshifted previous values of the predicted variable y, were used as input to predict the current value of y. The shortest time lag of x with a significant improvement over using only the values of y was deemed the effective lag length. Bonferroni correction to the level of significance (0.05) was applied to tests of each x-y pair.

Figure 3 shows the time lags of variables predicting the proportion of mask wearing and the difference in daily distancing at the county level, which provide a finer-grain view of the relationships than at the state level. For mask wearing, the most common lag length for R_t is 7 days, which could be partly attributed to the incubation period of COVID-19. Other more observable epidemiological indicators, such as COVID-like and influenza-like illnesses had a slightly shorter time lag. For daily distancing difference, mass media indicators had a wider influence in the shorter time frame, whereas R_t and the illness indicators had more uniform influence over the longer time lags. Note that not all time series variables in the data set appeared in the figure, because some of them, notably the social media variables, did not exhibit a significant lag with respect to the target variables (mask wearing and social distancing) according to Granger analysis.

Where applicable the variables were shifted by their respective effective lag lengths to account for the time lags between the variable and the behavior variable of interest. All subsequent analyses were performed using the time-shifted data.

The most relevant features were selected using forward stepwise regression on the shifted data. Starting with an empty set of predictors, at each step a linear regression model for predicting the target variable was constructed using the predictors selected so far and one of the unselected predictors. The models were evaluated using 10-fold cross validation, where the data set was randomly divided into 10 folds, and for each trial 9 folds were used for training and the remaining 1 fold was used for testing. The unselected variable that gave rise to the model with the highest average R² was added to the set of selected variables.

Figure 4 shows the R² of successive models with increasing number of predictors selected. The optimal cutoff was obtained by visually inspecting the graph and finding the “shoulder” where adding more variables did not result in much more improvement in the R² score. For mask wearing, the optimal cutoff was set at 6 variables. For daily distancing difference, the optimal cutoff was less clear. We selected two potential cutoff points, at 3 and 9 variables, for further exploration. The variables selected in each case are shown in Table 1.

For both mask wearing and daily distance difference, the predictor with most weight was a variable representing the political leaning towards Trump. For mask wearing, the additional variables spanned a variety of dimensions, including personality, epidemiology, media and demographics. For daily distance difference, the 3-variable set included additionally COVID-like illness and Population, whereas the 9-variable set was more heavily focused on the COVID-like and influenza-like illness indicators. As shown in Supplementary Figure S1 in the Supplementary Materials, these indicators were highly correlated. From Figure 4b we observed that the additional 6 variables gave rise to a noticeable though modest increase in R² over the 3-variable model. Overall, the best R² achieved by the mask wearing models were considerably larger than that of the daily distance models, indicating that the mask wearing model was a better fit than the daily distance model.

Our R-PVAs utilize a decision-making modeling approach known as Instance-Based Learning (IBL) (17, 18). Such models are particularly effective in dynamic and uncertain environments. The R-PVA models implementing IBL were developed in pyACTUP 2.2.3, which simulates the declarative memory module of ACT-R, including the storage and retrieval of experiences and a memory blending process. IBL leverages the ACT-R framework to simulate human decision-making based on past experiences stored in declarative memory. This memory stores multi-attribute past instances of decisions with subsets of attributes representing contexts (situational features), actions (decision; behavior choice), outcomes (resulting changes), and rewards. When a decision needs to be made, IBL uses a memory retrieval process called blended retrieval (or blending) to compute a choice or behavior based on decision context. Rather than retrieving a single best-matching past memory instance, blending computes a weighted average of multiple stored memory instances. The contribution of each instance to the blended result is determined by the activation, which is based on frequency and recency of the instance, and the similarity of the instance contexts to the decision context. Blending makes context-sensitive decisions that integrate past experiences and allows for generalization and adaptation in dynamic environments.

The mathematical foundation of IBL is rooted in Statistical Learning Theory (19), where it functions as a linear smoother—a non-parametric, instance-based learning function approximator. Consequently, IBL supports various learning paradigms, including Supervised Learning for regression and classification, as well as Reinforcement Learning for utility-based habitual behaviors. IBL continuously retains its “training data” within its memory repository, allowing it to adapt dynamically to new situations. To handle large datasets efficiently, IBL computations can be vectorized and parallelized, incorporating techniques such as approximate k-nearest neighbors for scalability. This combination of scalability, adaptability, and a cognitively plausible approach to decision-making through approximate expectations makes IBL particularly well-suited for modeling evolving human behavior in large-scale multi-agent systems.

ACT-R has subsymbolic mechanisms that determine the dynamics of the R-PVAs. Equations 1–3, below, define how the levels of activation of chunks in memory determine the probabilities of chunk retrieval.

Blended retrieval determines the value V that minimizes the sum of squared dissimilarities with the answer V_i proposed by each chunk, weighted by the probability, P_i, of retrieval of value V_i:(1)V=argminV∑iPi(1−Sim(V,Vi))2

The probability of retrieval is(2)Pi=eAi/s∑jeAj/sWhere the activation, A_i, is(3)Ai=Bi+∑fMPfSim(f,V)+εiand s and ɛ are noise factors, B is a base-level activation, and MP. MP_f is a mismatch penalty representing the dissimilarity between the representation of two values. Equation 4 defines how activation levels are increased by repeated experiences, or decay with time.(4)Bi=ln(∑j=1ntj−d)+βiwhere t_j is the time since the j^th storage or retrieval trial of chunk i, n is the number of occurrences, 0 < d < 1 is a decay parameter, and β_i is a constant offset. The parameters were set as follows: MP = 30.0, d = 0.5, and β = 0.0. The mismatch penalty MP was set to a value that allowed a relatively broad range of inexact matches to contribute to the activation level, weighted by the magnitude of the mismatch. The decay parameter d was set to the mid-point of the range, allowing for a moderate rate of decay.

Chunks generally can be represented as an unordered feature-value list of the form

{<feature₁: value₁>, <feature₂: value₂>, …, <feature_n: value_n>}

To make a prediction of the prototypical value of a feature given a (partial) set of predictor values, chunks that are similar to this chunk are retrieved and blended, weighted by a similarity function. Our R-PVA models used the similarity function(5)Sim(x,y)=1/(1+exp(y−x))2,withy>x.Thus, the mask wearing proportion of the next day, t + 1, was obtained by the blended retrieval of mask_wearing, given the values of the predictor variables on day t + 1:

{<smoothed_whh_cmnty_cli: v_t + 1>, <PctTrump_State_2016: v₁>,…, <TOM: v₆>},

In general, behavior change involves the gradual development of habits that are woven into the fabric of everyday life (20). A long history of psychological and neuroscience research on learning, decision-making, and behavior has led to the idea that there is a dichotomy between decisions and behaviors that result from effortful goal-striving or reflective processes vs. behaviors that are more automatic, less effortful habits (21–24). The execution of novel behaviors typically involves an initial goal-striving phase that requires cognitive effort in deliberation and behavioral control in the relevant environmental contexts. Then, repeated practice produces habit formation and strengthening processes that associate those specific behaviors to cues in the environment. Habitual behaviors come to be executed without effortful goal-striving cognitive processes. Thus, generally, novel behaviors are more difficult than familiar everyday behaviors because they require goal-striving.

The construct of self-efficacy in Social Cognitive Theory is a person's judgment about their ability to control events or their own behaviors. In general, an individual's belief in their ability to perform a behavior is necessary to attain the behavior. Self-efficacy has been defined in terms of an underlying cognitive learning process that may come from experienced mastery of the behavior or vicarious observation of similar others performing the behavior (14, 25). In our models, we assume mask-wearing is for most people (e.g., those not in health care) a novel behavior that initially requires effortful goal-striving processes, and that self-efficacy processes and repeated practice gradually make the behavior more likely to be executed.

For R-PVAs that model self-efficacy, two additional features are included:

{<predictor₁: value₁>, …, <predictor_n: value_n>, <difficulty: δ>, <effort: e>, <mask_value: m>},

In the current model, a “norm” is simply an initial behavioral bias in a subpopulation. During the initialization phase, a norm was established as the baseline for the model, using only the static variables and the values of the target (time series) variable from the initial 10 days of the data. For each state, the blended norm x of the target variable was obtained and 10 new chunks in proportion to this value x were learned, such that, keeping the other variable values as given, a proportion of x of the new chunks had a value representing the presence of the behavior (e.g., wearing a mask), and the rest of the 10 chunks had a value representing the absence of the behavior (e.g., not wearing a mask).

The self-efficacy mechanism in R-PVA involves two parameters, the initial effort needed to perform a behavior and the boost factor when the behavior is reinforced (without loss of generality, relative to a fixed level of difficulty and a threshold value that is a function of the other parameters). We performed a grid search over the parameter space to identify the optimal parameter setting, with the difficulty level set at 2.0 and using only data from California and Wyoming during the first wave of COVID-19 (up to 2020/06/30). The R-PVA models constructed were evaluated using n-fold rolling origin cross-validation (see Section 6.1).

Figure 5 shows the RMSE values of models for mask wearing and daily distance difference respectively. For mask wearing, the best-performing parameter values were: boost factor = 0.02, effort = 1.0. For daily distance difference, the best-performing parameter values were: boost factor = 0.01, effort = 1.8. We noted that for daily distance difference, the optimal parameter setting was at the boundary of reasonable parameter space, which might indicate that the self-efficacy mechanism was not of high utility. We will discuss this further below.

Norms obtained as described in Section 5.3 provided a starting point for the models. At each time point, ten memory chunks were added in proportion to the presence or absence of the target behavior as predicted by a R-PVA model. The model then simulated the behavior for the next time step using the present external variables and a blending of the past history of the internal memory chunks up to that point, subject to a decay function. The behavioral characteristics and the new memory chunks at the next time point were then set according to the new blended value.

For each of the two behavioral variables, mask wearing and daily distance difference, we constructed R-PVA models using the variables selected from the stepwise regression procedure and the parameter values obtained from the grid search. Several models of increasing complexity were evaluated, using a subset of the selected variables: •

Model 1: only the most significant time-series variable

Because of the sequential nature of the data, regular cross validation, with random assignment to folds, was not appropriate, as this would enable the prediction of a data point using future data points that should not have occurred yet. We instead analyze the models using n-fold rolling origin cross-validation (26), n being the length of the time series minus 1. For the i-th fold, the data from the initial time sequence <t_0, t_1,…t_i−1> was used for training, and the data at time point t_i was used for testing. Each successive training data set was a longer time series and included the previous training set. n-fold rolling origin cross-validation was performed using data in the date range 2020/04/24–2021/03/31. R² and RMSE scores were obtained for each of the models.

The evaluation results are shown in Tables 2, 3. For both mask wearing and daily distance difference, the models performed better as we included more of the selected variables (from R-PVA-*-1 to R-PVA-*-3), in particular going from the first model to the second model, suggesting that political leaning (PctTrump State 2016 and trump_approval_feb2020) had a substantial effect on the target behaviors.

The models of the two behaviors diverged with the inclusion of the self-efficacy mechanism. For mask wearing, self-efficacy helped improve the model, whereas for daily distance difference, the models with self-efficacy performed slightly worse than their corresponding models without self-efficacy [e.g., R-PVA-distance-3(a) vs. R-PVA-distance-4(a)]. We will discuss this in more detail in Section 7.

In addition, in the 9-variable setting for daily distance difference, both models without and with self-efficacy were worse than their 3-variable counterparts, suggesting that the 9-variable set might be overfitting the model. There were other indications that this might be the case. For instance, we noted earlier, when performing feature selection using stepwise regression (Figure 4b), that the addition of the 6 variables beyond the first 3 provided only a modest increase in R². Another consideration is that perhaps daily distancing was inherently harder to predict than mask wearing, as suggested by the lower maximum R² scores achieved for the two behaviors (Figures 4a, b). Note however, that R-PVA and regression models are based on different modelling principles and the results from stepwise regression, although suggestive, are not directly extensible to R-PVAs.

We noted earlier that the variables that carried the most weight in the feature selection procedure were political variables (“PctTrump State 2016” and “trump_approval_feb2020” respectively), suggesting that political partisanship is a strong factor in the adherence to non-pharmaceutical interventions during COVID-19. Here we will delve into this aspect more closely by examining the states with the most extreme values for “PctTrump State 2016”. The five states with the highest “PctTrump State 2016” values were West Virginia, Wyoming, Oklahoma, North Dakota and Kentucky. The five states with the lowest “PctTrump State 2016” values were Maryland, New York, Hawaii, Vermont and California. We will refer to these two sets of states as HighTrump and LowTrump states respectively.

Figure 6 shows the observed and predicted values of the two behavioral variables in the HighTrump and LowTrump states over time. The predicted values were obtained from the best performing models of the two behaviors (R-PVA-mask-4 and R-PVA-distance-3). In general the predicted values tracked the observed values fairly closely. LowTrump states showed a higher proportion of mask wearing behavior and a larger decrease in distance travelled compared to the pre-COVID baseline. For mask wearing, the values in HighTrump states fluctuated but exhibited an overall upward trend over time, whereas the values in LowTrump states showed less fluctuation, possibly because the proportion of mask wearing was already fairly high in those states. For difference in daily distance, the values in both HighTrump and LowTrump states displayed more variations than those for mask wearing, perhaps another factor contributing to the lower performance of the social distancing models, when compared to the mask wearing models. For most states, the difference in daily distancing progressed in an S-shape, increasing steadily till around day 80 and then decreasing till around day 280 when it started to rise again, drastically in some cases. These inflection points roughly corresponded to the end of the first wave (and start of the second wave) of COVID-19 and the vaccine roll out at the end of the third wave. The HighTrump states had smaller decreases in distance travelled, but with slightly more variation than LowTrump states, and in some states there was even an increase in distance travelled compared to the pre-COVID baseline.