---
title             : "Studying behaviour change mechanisms under complexity"
shorttitle        : "Complexity in behaviour change mechanisms"
author: 
  - name          : "Matti T. J. Heino"
    affiliation   : "1"
    email         : "matti.tj.heino@gmail.com"
  - name          : "Keegan Knittle"
    affiliation   : "1"
    email         : "keegan.knittle@helsinki.fi"
  - name          : "Chris Noone"
    affiliation   : "2"
    email         : "chris.noone@nuigalway.ie"
  - name          : "Fred Hasselman"
    affiliation   : "3"
    email         : "f.hasselman@pwo.ru.nl"
  - name          : "Nelli Hankonen"
    affiliation   : "1"
    corresponding : yes    # Define only one corresponding author
    email         : "nelli.hankonen@helsinki.fi"
    address       : "Faculty of Social Sciences, University of Helsinki, PO Box 54, 00014 University of Helsinki, Finland"
affiliation:
  - id            : "1"
    institution   : "Faculty of Social Sciences, University of Helsinki, PO Box 54, 00014 University of Helsinki, Finland"
  - id            : "2"
    institution   : "School of Psychology, National University of Ireland, Galway, Ireland"
  - id            : "3"
    institution   : "Behavioural Science Institute, Radboud University Nijmegen, The Netherlands"
note              : "This is a preprint submitted for publication"
abstract: | 
  Understanding the mechanisms underlying the effects of behaviour change interventions is vital for accumulating valid scientific evidence, and useful to informing practice and policy-making across multiple domains. Traditional approaches to such evaluations have applied study designs and statistical models, which implicitly assume that change is linear, constant and caused by independent influences on behaviour (such as behaviour change techniques). This article illustrates limitations of these standard tools, and considers the benefits of adopting a complex adaptive systems approach to behaviour change research. It 1) outlines the complexity of behaviours and behaviour change interventions, 2) introduces readers to some key features of complex systems and how these relate to human behaviour change, and 3) provides suggestions for how researchers can better account for implications of complexity in analysing change mechanisms. We focus on three common features of complex systems (i.e. interconnectedness, non-ergodicity and non-linearity), and introduce Recurrence Analysis, a method for nonlinear time series analysis which is able to quantify complex dynamics. The supplemental website (https://git.io/Jffrm) provides exemplifying code and data for practical analysis applications. The complex adaptive systems approach can complement traditional investigations by opening up novel avenues for understanding and theorising about the dynamics of behaviour change.  

keywords: "complex systems, wellbeing, methodology, behaviour change"
linenumbers       : yes
figsintext        : yes
floatsintext      : yes
figurelist        : no
tablelist         : no
footnotelist      : no
lineno            : yes
mask              : no
csl               : behavioral-sciences.csl
bibliography      : ["complexity-behchange_3.bib"] # "r-and-packages.bib"
lang              : "en-GB"
documentclass     : "apa6"
classoption       : "man"  # man / doc / jou (& other comma separated options)
output: 
  papaja::apa6_pdf: default
  pdf_document    : 
    fig_caption: true
  # papaja::apa6_word:
    # reference_docx: behavsci.docx
header-includes:
  - \usepackage{booktabs}
  - \usepackage{float}
  - \usepackage{setspace}
  - \AtBeginEnvironment{tabular}{\singlespacing}
  - \AtBeginEnvironment{lltable}{\singlespacing}
  - \AtBeginEnvironment{tablenotes}{\doublespacing}
  - \captionsetup[table]{font={stretch=1.5}}
  - \captionsetup[figure]{font={stretch=1.5}}
---

```{r load-data, include = FALSE}

library(knitr)
library(tidyverse)
library(bookdown)

knitr::opts_chunk$set(echo = FALSE, 
               warning = FALSE,
               message = FALSE,
               error = TRUE,
               cache = FALSE, 
               collapse = TRUE,
               dpi = 300
               )



# rmarkdown::render(input = "_complexity-manuscript.Rmd", output_format = "papaja::apa6_word", encoding = "UTF-8")
# citr::tidy_bib_file(rmd_file = "_complexity-manuscript.Rmd", messy_bibliography = "complexity-behchange_3.bib")
```

\newpage

# Introduction

In order to understand why behavioural interventions often fail to produce sustainable effects [@kwasnickaTheoreticalExplanationsMaintenance2016], especially when transferred from one context to another, a core interest of behaviour change science is to improve our understanding of mechanisms of behaviour change. Behavioural theories identify hundreds of potential "determinants" of behaviour, that is, factors that potentially influence the behaviour of interest, constituting the mechanisms by which behaviour change techniques might influence behaviour [@careyBehaviorChangeTechniques2019]. These range from social cognitions such as self-efficacy and attitudes, to biological factors, and certain elements of the social and built environments in which behaviours take place [@michieABCBehaviourChange2014]. When studied using typical factorial designs and linear statistical models, the relationships between causal precedents and behaviour change are assumed additive, constant and linear (i.e. the outputs are proportional to the inputs). However, it is our position that this offers behaviour change researchers and the general public an inaccurate or at least imprecise understanding of behaviour change. We should consider the relevant factors as complex, potentially non-linear, and dynamic. 

The evaluation of behaviour change interventions often involves randomly assigning participants to receive an intervention of interest or a specific comparator and measuring subjective and objective indicators of behaviour [@matthews2020EvaluationBehaviorChange]. Usually, these measurements occur immediately before and after the delivery of the intervention, though sometimes additional follow-up measurements may take place weeks or months later. This is the classic Randomised Controlled Trial design and the data produced are most often analysed using statistical techniques that are specific cases of the General Linear Model. In this paper, we refer to this as the “conventional approach”. If the interest is only in assessing whether the treatment overall was more effective, on average, in the intervention group than the control group, comparing averages in randomised controlled trials can be purposeful and acceptable (i.e. answering questions such as "Does the intervention have an effect on the target behaviour?", "Do cohorts differ from each other?"). However, using so few measurement points to study behaviour change mechanisms ("How do intervention participants change?") may present problems. 

Currently, mechanisms of change in behaviour change interventions are typically studied using mediation analysis [see e.g. @haggerKnownKnownsKnown2020], where one variable’s (X) impact on another (Y) is modelled to pass through a third variable (M). In its classical form, one expects the path X-Y to go near zero when adding M to the model. If this is observed, the researcher concludes there is evidence for mediation. For example, within an RCT, researchers can collect not only outcome data (e.g. physical activity behaviour, Y), but so called process variables M (e.g. autonomous motivation for physical activity, self-regulation skills), reflecting e.g. psychological mechanisms that are hypothesised to be the explanation for what makes the intervention (X) effective in changing behaviour (case example from members of the current author group: [@hankonenWhatExplainsSocioeconomic2017]). 

Inferring mechanisms from contemporary mediation analyses is problematic on various grounds [@bullockYesWhatMechanism2010; @fiedlerWhatMediationAnalysis2011; @greenEnoughAlreadyBlack2010; @klineMediationMyth2015; @roeWhatWrongMediators2012; @keele2015CausalMediationAnalysis], but of particular importance to human behaviour change, the results depend on four key assumptions, discussed by [@hofmann2020LinearMediationDynamic]. We can only consider the results of such mediation analyses to be accurate when we assume that: 1) The number of variables involved is small, and dynamics can be meaningfully assessed with only a few time points, 2) The process of change is the same for all individuals, e.g. follows the same sequence, 3) The dynamic between variables is linear, additive, and does not change in time, and 4) The included variables are not entangled with the context, omitted variables, or each other, in bi-directional recursive relationships. Researchers can of course include more variables (leading to new issues, e.g. mixing up mediators, confounders and colliders [@rohrerThinkingClearlyCorrelations2018] or lowering the already worrisome statistical power [@loeysCautionaryNotePower2015; @schoemannDeterminingPowerSample2017]), try to incorporate non-linear effects [@knafl2017IncorporatingNonlinearityMediation], and add more time points (for caveats regarding latent variable modeling, see e.g. [@hamakerCritiqueCrosslaggedPanel2015; @pmlr-v6-dawid10a]). But limiting the notion of a mechanism to a (multiple) mediation/moderation problem limits our understanding of how changes occur over time [@bechtel2013ThinkingDynamicallyBiological; @golonka2019EcologicalMechanismsCognitive]. 

During the first two decades of the 21st century, behaviour change researchers have started embracing designs with an increased focus on temporal processes [@bolgerIntensiveLongitudinalMethods2013; @reichert2020AmbulatoryAssessmentPhysical], extending the traditional approach. Recently, alternative solutions stemming from complex systems science [@siegenfeldIntroductionComplexSystems2019; @newmanComplexSystemsSurvey2011] have become increasingly accessible and helpful in tackling problems of understanding change processes. We will explore these ideas and how they help us surpass traditional assumptions. In what follows, this paper will 1) outline the complexity of behaviours and behaviour change interventions, 2) introduce readers to some key features of complex systems and how these can be applied to human behaviour, and 3) provide concrete suggestions for how researchers can better account for the implications of complexity in analysing behaviour change mechanisms. 

## What are complex systems?

A system is "a delineated part of the universe which is distinguished from the rest by an imaginary boundary" [@bar-yamConceptsSystem2018], although other definitions exist [see @wrightThinkingSystemsPrimer2009 for a primer]. Many things---an airplane, a car, a robot, a central nervous system, a school, a community, a society---can be conceptualised as systems. This paper focuses on individual people, which are complex systems. Complex systems can be characterised as webs of many _interdependent self-organising_ parts that operate without central control, whose interactions give rise to _emergent_ properties and behaviours, irreducible to a sum of its parts [@mitchellComplexityGuidedTour2009]. The future behaviour of such a complex system strongly depends on its _unique history of interactions_, that is, past experience. Additionally, the system adapts to each its environment and actors therein, _coevolving_ with each other to create macro-level behaviour, which is difficult to predict and usually not changeable in a stepwise engineering sense [@brandTailoringHealthyWorkplace2015]. These characteristics (in italics) distinguish complex systems from those which are just complicated. Highly complicated processes or systems (e.g. an airplane), unlike complex ones (e.g. an organism) cannot, for example, self-organise to function adaptively when a part is removed [@ricklesSimpleGuideChaos2007]. Guides to basic terminology of complexity for scientists working with health behaviours can be found in [@ricklesSimpleGuideChaos2007] as well as table 1 of [@brandTailoringHealthyWorkplace2015]. The promise of complex adaptive systems approaches in health behaviour change research has been previously discussed by e.g. [@resnicowChaoticViewBehavior2006], but over a decade ago, not many empirical solutions were easily accessible to investigators in this field. Recently, methods to study complex adaptive systems in health behaviour change have been presented @gomersallComplexAdaptiveSystems2018, with a focus on simulation and qualitative methods. This article discusses novel quantitative solutions, which have recently become available, to investigate behavior change phenomena with a complex systems lens.

## The relevance of complexity for behaviour change

To paint a picture of just how complex the behavioural world is, take the case of physical activity as an example behaviour. Already three and a half decades ago, more than 30 influences on (or "determinants of") this behaviour were being considered, along with calls for better understanding of their dynamics, interactions, and the time scales over which these develop [@dishmanDeterminantsPhysicalActivity1985]. While any influence (e.g. intention, attitude) could have a direct relationship with physical activity, some rely on interactions with other influences to affect behaviour [e.g. preventive behaviours being dependent on fear only in the presence of sufficient efficacy beliefs; @kokIgnoringTheoryMisinterpreting2018; @petersConsensusFearAppeals2018]. Furthermore, these interactions may be moderated by additional factors, and by other variables which themselves have no direct relationship with physical activity, with synergistic and opposing effects which may themselves depend on whether some threshold is exceeded. The extent to which all known (and unknown) influences on physical activity interact with one another presents a map of practically infinite, intertwined ‘routes’ to initiating and maintaining the activity. 

The role of time brings added complexity to this behavioural world, as dynamic patterns change over time and at varying frequencies [@scholz2019ItTimeThink; @bolger2019HeterogeneityTemporalProcesses]. For example, fluctuations in physical activity clearly occur within a day, as most individuals are (at least in the absence of highly sedentary working conditions and considerable somnambulism) more active while awake than while asleep. Fluctuation also occurs over the course of a week, as activity levels tend to be higher on weekdays than on weekends [@matthewsSourcesVarianceDaily2002]; over the course of months, as activity levels are higher in warmer seasons and lower in colder ones [@cepedaSeasonalityPhysicalActivity2018]; and over the course of years, as activity levels tend to decline with age [@dumithPhysicalActivityChange2011]. How determinants---which are postulated to comprise the mechanisms underlying changes in behaviour---fluctuate and interact with the fluctuations in behaviour, is largely unknown.

Human behaviour is complex, and while we have formulated theoretical constructs to be as amenable as possible to linear methods of analysis, this may obscure important characteristics of behaviour change. Why are linear models inappropriate for many of our research questions in the behavioural sciences? First, with many nonlinear interactions across time scales, mechanistic causality (including mediation and moderation) becomes suspect or intractable [@richardsonInteractionDominantDynamicsTimescale2017; @ricklesCausalityComplexInterventions2009; @pmlr-v6-dawid10a]. Second, traditional statistical analyses start from the simplification that everything is independent from everything else, whereas in actuality, nearly everything eventually depends on everything else, contributing to what Paul E. Meehl [-@meehlWhySummariesResearch1990, p. 204] seminally coined as "the crud factor". Jacob Cohen, the developer of power analysis, similarly exclaimed that [in the absence of randomisation], the nil hypothesis of no effect is always a priori false [@cohenEarth051994] -- these well-known ideas demonstrate violations of the classical assumptions regarding independence and interference [@finkSocialDeterminantsPopulation2016; @wallotInteractionDominantCausationMind2017]. In the same vein, forecasting in complex systems is notoriously difficult [@makridakisDecisionMakingPlanning2009; @makridakisForecastingSocialSettings2019; @taleb2020StatisticalConsequencesFat], making hypothesis testing---which is, after all, the test of a prediction---in intervention evaluation a curious challenge. Complexity science, which starts from the assumption that everything is intertwined, can provide us with new hypotheses which respect the complexity of the phenomena under study [@heinoCommentaryReproducibilityPsychological2017]. This is necessary, because a conventional linear analysis will only give results that are correct given the assumption that the components in the model are independent, with additive effects that can be decomposed and attributed to their causes (e.g. beta coefficients in multiple regression). If, on the other hand, these "component-dominant" dynamics are not driving the system, but instead the effects are intertwined, overlapping and inseparable [as proposed in the behaviour change context by @petersPragmaticNihilismHow2017], and thus the dynamics are "interaction-dominant", then replication and generalisation issues for results stemming from the linear analysis are almost inevitable [@wallotInteractionDominantCausationMind2017]. Instead, intensive longitudinal methods are necessary to monitor how processes unfold; this information can then be used to dynamically tune interventions in real time, making success less dependent on having a correct program theory at the outset [@heklerWhyWeNeed2019; @schiepekRealtimeMonitoringPsychotherapeutic2016]. Made possible by N-of-1 methodologies [@kwasnickaChallengesSolutionsNof12019], this goal has been recently pursued by using e.g. control systems engineering approaches [@heklerTutorialUsingControl2018] and just-in-time interventions [@hardeman2019SystematicReviewJustintime].

Although e.g. behaviour change maintenance has been theorised at length [@kwasnickaTheoreticalExplanationsMaintenance2016], theories outlining complex systems principles (which largely avoid the aforementioned issues) have been missing from this work. From the viewpoint of complex systems science, the effects of behaviour change interventions can be considered as shocks to the system in which they take place -- the aim of the shock is to alter the system’s status, pushing against existing forces to affect change [@haweTheorisingInterventionsEvents2009; @olthofDestabilizationSelfratingsPsychotherapeutic2019]. This is akin to attempts to work against gravity, which pulls a ball in a valley (a relatively stable state, also known as an attractor; see Figure \@ref(fig:ballvalley)) to the bottom of it [@schefferQuantifyingResilienceHumans2018, p. 11885; @ashwin2012TippingPointsOpen]. Taking the analogy further, pushing the ball outside of the valley may lead it to roll off a peak, ending up in a deeper valley (i.e. less unstable, more deep-rooted state) than from where it started. A complex systems perspective implies, that even in the event of a successful intervention, stabilizing a system in a more functional state may require at least as many resources as the initial change itself [@bar-yamMakingThingsWork2004, p. 211].  In general, while complex systems may often be impossible to control precisely, they can be stewarded approximately, while allowing for variability stemming from self-organisation to flourish instead of trying to iron it out [@navarroHealthyVariabilityOrganizational2015; @talebBlackSwanCairo2011]. The necessity of complex systems approach is increasingly recognized; for example, it is highlighted in the UK Medical Research Council’s recently updated guidance for development and evaluation of complex interventions [@skivingtonDevelopingEvaluatingComplex2018].

```{r ballvalley, fig.width = 8.27, fig.height = 11.69/6, fig.cap = "Evolution in attractor landscape: An intervention moulds a system, making it less stable, hence easier for the ball to move from current state (left) to another one (right). Alternatively, an intervention -- or random events -- can jolt the system over the ridge, i.e. a tipping point."}

# knitr::include_graphics("figures/ball_valley.png")

grid::grid.raster(png::readPNG("./figures/ball_valley.png"))

```

Having now undergone a brief conceptual introduction to complexity, we can describe behaviour change as _"a collection of contextualised processes that are nontrivially specific to each individual, and which form a complex interconnected system that is not restricted to linear dynamics"_ (see [@wrightPersonalizedModelsPsychopathology2020a], p. 4). We highlight three features of this definition:

1. _A complex interconnected system_: A multitude of variables and timescales which are interwoven, interdependent, and interacting.

2. _Contextualised processes, specific to each individual_: Individuals follow meaningfully different change trajectories that develop and change with time.

3. _Not restricted to linear dynamics_: Inputs are not necessarily proportional to outputs, and long periods of apparent stability can precede short periods of rapid change.

# Behaviour change mechanisms under complexity: Three key features

In the following three sections, we drill further down into these ideas. In the first, we introduce interconnectedness via interaction-dominant dynamics, which flow from point 1 above; second, we present how idiosyncratic, non-stationary change trajectories lead to non-ergodicity, a technical term for point 2; third, we highlight that the flexibility of complex systems leads to ubiquitous nonlinear dynamics as alluded to in point 3. Table 1 provides an overview of these ideas, which are elaborated further in the subsequent sections. 

```{r summary-table, results = "asis"}
text_tbl <- dplyr::tribble(
  ~" ", ~"Interconnectedness", ~"Non-ergodicity", ~"Non-linear dynamics",
  "Description",
    "The structure of a system---how it is organised and the relationships between its component parts---can matter more than the component parts themselves. This includes interconnectedness of different variables such as attitudes or perceived norms, as well as that of their temporal dependence; dynamic dependencies of complex systems are not restricted to one or a few previous time points. @richardsonInteractionDominantDynamicsTimescale2017 @wallotInteractionDominantCausationMind2017 @vanrooijFractalApproachDynamic2013",
    "Psychological processes are non-stationary and heterogeneous, hence non-ergodic (group-level measurements do not correspond to those of individuals in time). This means within-individual processes cannot, in general, be inferred from between-individual data. The lack of group-to-individual generalisability implies a threat to validity of results in many areas of science. @fisherLackGrouptoindividualGeneralizability2018 @molenaarImplicationsClassicalErgodic2008 @2019TimeMoveAverage",
  "In a linear progression of a phenomenon, the whole is exactly the sum of its parts: You can calculate how much each influencer of behaviour changes, and add them together to get the total effect. Non-linearity occurs when a systems’s inputs are disproportionate to its outputs. For example, an effect might be imperceptible for a long time, then explode (as in exponential growth), or suddenly switch states upon reaching a threshold. @helmichSuddenGainsDaytoday2020 @kelty-stephen2017MultifractalityMonoFractality @olthof2020ComplexityPsychologicalSelfratings",
  "Main lesson", 
    "Dynamic, intertwined processes do not exist in a vacuum, but are always co-dependent and cannot be partialed out into variance components without losing essential information on how the system as a whole operates.",
    "Drawing individual-level inferences from group-level data (the ecological fallacy) leads to misleading or incorrect inferences regarding individual behaviour. A statistical relationship in the population may not hold for any of the individuals.",
    "Viewing the world solely from the lens of linear phenomena and relationships, leads to missed opportunities and misunderstood impacts of interventions.",
  "Recommendations for the research community",
    "Move from traditional regression-based approaches, which are inspired by component-dominant, additive dynamics (that the effects, or beta coefficients, of each variable can be summed together), linear approximations and Gaussian distributions, to methods able to cope with multiplicative effects and heavy-tailed distributions.",
    "Move from large-sample research with many variables and many people but few time points (one model per sample), to N-of-1 and intensive longitudinal time series designs, with usually fewer people and variables, but more data per variable (one model per individual).",
    "Move from linear approximations with the illusion of predictability, to methods that can accommodate non-linear patterns and disproportionate influences.")

# knitr::kable(text_tbl)

papaja::apa_table(text_tbl, caption = "Three common features of complex systems, with recommendations for behaviour change research.")
```

## Interconnectedness

When processes in complex systems are not independent, they are said to be coupled. Coupling can be unidirectional (where, for example, physical activity increases muscle mass but not the other way around), or bidirectional, where the elements of a system (e.g. good performance and rewards) simultaneously reinforce or suppress each other as time progresses, demonstrating a type of circular causality. Dynamics in living systems tend to be dominated by synergies ("interaction-dominant causation") instead of their component parts ["component-dominant causation"; @bakSelforganizedCriticalityExplanation1987; @richardsonInteractionDominantDynamicsTimescale2017; @wallotInteractionDominantCausationMind2017]. Many psychological and behaviour change theories seem to at least implicitly assume the presence of reciprocal causation and intertwined processes (e.g. [@banduraSocialFoundationsThought1986], p. 6), but empirical testing of such processes has to date been limited. 

As mentioned earlier, within the conventional approach to behaviour change intervention evaluation, researchers commonly employ mediation analyses to examine mechanisms. However, given its assumptions, the clean _independent variable_ $\rightarrow$ _mediator_ $\rightarrow$ _dependent variable_ type of path analysis can be misleading, when change is in fact driven by self-reinforcing, "autocatalytic" interactions, flowing in time. In component-dominant causation, effects follow causes in a billiard-ball fashion, and one variable can change without everything else changing. For example, a study developed with the component-dominant mindset could aim to find out how using a specific behaviour change technique, say goal setting, affects behaviour. On the other hand, variables of interest to behaviour change researchers are unlikely to change without affecting a large amount of other, related variables [@petersPragmaticNihilismHow2017], producing highly context-dependent effects [@craigTakingAccountContext2018]. This, too, implies that interaction-dominant causation is a more plausible framework for the behaviour change domain, wherein effects emerge (and are conditional upon) the system’s holistic multivariate dynamics, with everything potentially taking place simultaneously in a circularly causal manner. Interaction-dominant dynamics are also characterised by heavy-tailed distributions [@vanrooijFractalApproachDynamic2013], such as the log-normal distribution [@trafimowMeansStandardDeviations2018], which are common in psychological data [@cainUnivariateMultivariateSkewness2017; @bonoNonnormalDistributionsCommonly2017], as well as the presence of long-range temporal correlations and power-law scaling [@vanordenHumanCognitionScaling2005; @wijnantsReviewTheoreticalPerspectives2014a; @olthof2020ComplexityPsychologicalSelfratings]. Importantly, interplay happens not just between variables, but also their temporal dynamics: Processes taking place on fast timescales (e.g. lack of physical activity) modulate slow-timescale processes (e.g. development of obesity, lower energy levels), which feed back and affect the fast-timescale processes [@richardsonInteractionDominantDynamicsTimescale2017]. 

One way of looking at mutually interacting processes with reciprocal causality is to consider the system as a network. Network science is a well-established field with applications ranging from physiology to the organisation of cities [@barabasiNetworkScience2016], and health [@centolaTruthBehavioralChange2018; @zhangSocialNetworksHealth2019]. An illustrative example comes from the study of depression, where the traditional latent variable thinking assumes that a latent factor—depression—causes the symptoms. On the contrary, a network science perspective leads to an alternative view, where the network of mutually interacting symptoms constitutes the phenomenon [@borsboomNetworkTheoryMental2017; @cramerMajorDepressionComplex2016]. This approach has provided new avenues into understanding and treating depression, such as locating the symptoms which are most relevant to the activation of the network (i.e. the emergence of depression), or considering how intervening on specific symptoms might affect the system, given all dampening and reinforcing pairwise relationships between symptoms.

Although the network theory of mental disorders [@borsboomNetworkTheoryMental2017] aligns with and stems from complexity science, the psychological network models usually associated with the approach [for applications in health psychology, see @heinoVisualisationNetworkAnalysis2019a; @mkhitaryanNetworkApproachHealth2019] rely on many assumptions stemming from their grounding in multiple regression; including multivariate normality (i.e. linearity) and stationarity [for a comprehensive treatment, see @epskampGaussianGraphicalModel2018], as well as being very different from their physical counterparts with properties such as nonlinear scaling and space-filling [@westHomeostasisGaussStatistics2010; @westScaleUniversalLaws2017]. Still, the conceptual frameworks such models represent---coupled processes interacting in a system, instead of "root causes" [@bringmannDonBlameModel2018a]---ought to be the primary ontology considered by behaviour change researchers. In the later section on empirical solutions, we present a recurrence-based network modelling approach to consider these coupled processes  [@hasselmanStudyingComplexAdaptive2020].

## Non-ergodicity

To be useful to individuals, processes postulated by psychology ought to work on the individual level [@johnstonUsefulTheoriesShould2013]. Whether group-level variation is informative of individual-level dynamics, depends on a condition known as ergodicity, which has the following properties: "Only if the ensemble of time-dependent trajectories in behavior space obeys two rigorous conditions will an analysis of interindividual variation yield the same results as an analysis of intraindividual variation [...] First, the trajectory of each subject in the ensemble has to obey exactly the same dynamical laws (homogeneity of the ensemble). Second, each trajectory should have constant statistical characteristics in time (stationarity, i.e., constant mean level and serial dependencies)" ([@molenaar2008ConsequencesErgodicTheorems]; see also [@molenaarImplicationsClassicalErgodic2008]).

In other words, this would mean that in a 100x100 spreadsheet, where participants are rows and measurement occasions are columns, calculating an average of values within one column ("ensemble average"), would give the same result as calculating the same statistic from one row ("time average"). For example, in an ergodic process, the mean and standard deviation of each person’s daily minutes of physical activity over a 100-day period would be the same as the mean and standard deviation of 100 people’s daily physical activity minutes measured once. Or, observing that 20% of a given population are smokers, would mean that everyone is a smoker for 20% of their lives. In terms of coupled processes, the correlation between physical activity and intention would be the same in the population measured once, as it is for one person over time. 

Going back to the two "rigorous conditions", the condition of homogeneity almost by definition rules out the behaviour change researcher’s interests, as we are interested in how people (can) change, and it is quite clear that people do not all follow the same behaviour change processes. Indeed, it would seem preposterous to suggest that, for example, self-regulation is a constant process during a individual’s life span. Although the mathematical proof for the non-equivalence of inter-individual and intra-individual data structures was published over a decade ago [@molenaarManifestoPsychologyIdiographic2004], only recently has serious research attempted to quantify the threat stemming from lack of group-to-individual generalisability [@fisherLackGrouptoindividualGeneralizability2018]. This preliminary work indicates that even if we could work with "generalisable" ideal random samples from well-defined populations, we would still be committing the ecological fallacy (i.e. drawing individual-level inferences from group-level data) if we wanted to apply our knowledge to individuals. 

The second condition, that the statistical properties of these processes must not change over time, is generally referred to as stationarity. In the context of physical activity, the extent to which physical activity is influenced by other factors, is likely to change over time. For example, the effect of discomfort on physical activity is likely to change in a non-linear manner over time, as fitness and tolerance of discomfort fluctuate not only because of randomness, but as core features of the phenomena itself [@hallTemporalSelfregulationTheory2015]. However, the tools most often used in research for thinking about and analysing behaviour change, such as linear regression, do not account for these kinds of temporal dynamics. This is because temporal cognitive change fundamentally violates the assumption of stationarity, as exemplified next.

For the processes underlying physical activity outlined above to be considered stationary, the average level of discomfort must remain stable across time for all individuals. In addition, the sequential dependence between repeated measures must be stable [i.e. the variance must be constant and the sequential correlations must only be influenced by how far away in time two data points are; @molenaarNewPersonspecificParadigm2009]. In terms of the relationships between variables, the assumption of stationarity requires that the causal structure which leads to a particular outcome is unchanging across time [@coleTestingMediationalModels2003]. Examining behaviour change usually involves an attempt to change the causal structure underlying a behaviour (e.g. after learning to make coping plans to tackle barriers to physical activity, the causal relationship from perceiving a barrier to subsequently deviating from one’s plan to be active, ought to be diminished), and generally means that either a decrease or increase in a particular behaviour is expected as learning and development progress. Stationary data is therefore rare in behaviour change research. This lack of stationarity has however rarely been acknowledged or (statistically) accounted for in empirical studies evaluating behavioural processes. The result is analogous to the ecological fallacy of taking a population-level mean and extrapolating to individual-level attributes; an average over an individual’s time series describes that individual better than the population-level snapshot, but still might not applicable to any particular time period. As a simple example, think 100 days of data in which a linear dependence relationship is strongly positive for the first 50, and strongly negative the other; you might observe the average correlation over the whole time series to be zero.

Figure \@ref(fig:tv-var) illustrates non-stationarity in the case of work motivation, a key feature of occupational health psychology. Data is from one participant in an observational study of motivation self-management (Heino et al., in prep). We can observe that the relationships vary drastically, as the study progresses.

```{r tv-var-prep, include = FALSE, eval = FALSE}

library(tidyverse)
library(mgm)

sdtVarData <- readr::read_rds("./data/EMA_data_Moti_P10.csv")
username <- "Moti_P10"

sdtVarData <- sdtVarData %>% 
  dplyr::mutate_at(vars(anxiety_guilt, autonomy, competence, 
                        for_others, importance, interest, own_question, pleasure, 
                        productivity_work, satisfaction_work, sleep_quality,
                        relatedness, required, time_since_break), 
                   funs(as.numeric)) %>% 
  # dplyr::mutate(`interesting_pleasant` = (interest + pleasure) / 2) %>%
  dplyr::select(# autonomy, competence, relatedness, 
                interest, pleasure, Importance = importance,
                `situation_requires` = required, `avoid_guilt_anxiety` = anxiety_guilt, 
                `another_wants` = for_others, 
                dateTime) %>% 
  dplyr::arrange(dateTime) %>% 
  dplyr::mutate(timediff = c(NA, diff(dateTime))) %>% 
  dplyr::filter(timediff > 15) %>% 
  dplyr::select(-timediff) %>% 
  na.omit(.)

df <- sdtVarData$dateTime
time_01 <- df - df[1]
time_01 <- as.numeric(time_01)
time_01 <- time_01 / max(time_01) # normalise time, not necessary for mgm, though

sdtVarData_values <- sdtVarData %>% 
  dplyr::select(-dateTime) %>% 
  dplyr::mutate_all(funs(as.numeric(.)))

bwSeq <- seq(0.01, 1, length = 10) # ten equally spaced values in [0:01; 1]

# Find optimal bandwidth
bw_object <- mgm::bwSelect(data = sdtVarData_values,
  type = rep("g", 6),
  level = rep(1, 6),
  bwSeq = bwSeq,
  bwFolds = 2,
  bwFoldsize = 20,
  modeltype = "mvar",
  lags = 1,
  scale = TRUE,
  timepoints = time_01)

readr::write_rds(bw_object, path  = paste0("RDS_files/", username, "_6nodes_bw_object.RDS"))
bw_object <- readr::read_rds(path  = paste0("RDS_files/", username, "_6nodes_bw_object.RDS"))

bandwidth <- bwSeq[which.min(bw_object$meanError)]

# data.frame(bwSeq, bw_object$meanError) %>% plot()

# --------- Fit time-varying VAR -----------------------------------

estpoints <- seq(from = 0, to = 1, length = 50)

set.seed(1)
tvvar_obj <- mgm::tvmvar(data = sdtVarData_values,
  type = rep("g", 6),
  level = rep(1, 6),
  lambdaSel = "CV",
  timepoints = time_01,
  estpoints = estpoints,
  bandwidth = bandwidth,
  lags = 1,
  saveData = TRUE,
  scale = TRUE)

readr::write_rds(tvvar_obj, path  = paste0("RDS_files/", username, "_6nodes_tvvar_obj.RDS"))
tvvar_obj <- readr::read_rds(path  = paste0("RDS_files/", username, "_6nodes_tvvar_obj.RDS"))

resample_obj <- mgm::resample(object = tvvar_obj,
  data = sdtVarData_values,
  nB = 50,
  blocks = 10,
  seeds = 1:50,
  quantiles = c(.05, .95))

readr::write_rds(resample_obj, path  = paste0("RDS_files/", username, "_6nodes_resample_obj.RDS"))
resample_obj <- readr::read_rds(path  = paste0("RDS_files/", username, "_6nodes_resample_obj.RDS"))

```

```{r tv-var, fig.width = 8.27, fig.height = 8.27/3, fig.cap = "Relationships between a single participant’s motivational variables varying in time (time-varying autoregressive model). Networks represent relationships between variables around the time points where 10% (panel A), 50% (B) and 90% (C) of the study had been completed. An arrow from one variable to the next means the former predicts the latter at the next time point; green for positive and red for negative correlation. If a stationary model was used, all periods would be collapsed to a single result, creating the impression that the relationships were homogeneous across the study period. Although this temporal variability can be due to e.g. changes in how the participant answers the questions (boredom, shifting perception of the items, etc.), or poor reliability of the measures, complexity theory would also guide us to expect that in very concrete reality, the direction and strength of relationships can shift over time and differ based on the state a person resides in. As an example, the relationships between motivational variables during behaviour change initiation phase, may differ from the relationships during the maintenance phase."}

username <- "Moti_P10"
tvvar_obj <- readr::read_rds(file = paste0("./RDS_files/", username, "_6nodes_tvvar_obj.RDS"))

# Maximum edge weight to scale edges to:

maxedge <- max(c(tvvar_obj$wadj[, , 1, 5], tvvar_obj$wadj[, , 1, 25], tvvar_obj$wadj[, , 1, 45]))

layout(matrix(
  c(1, 2, 3), 
  nrow = 1, ncol = 3, byrow = FALSE))

vertex_size <- 30
margin_size <- c(4.5, 5, 4.5, 5)
label_size <- 1.3
loop_rotation <- 0
title_size <- 1.5
  
# 10% of study completed
tp <- 5
qgraph::qgraph(t(tvvar_obj$wadj[, , 1, tp]), 
               layout = "circle",
               edge.color = t(tvvar_obj$edgecolor[, , 1, tp]),
               curveAll = FALSE, parallelEdge = TRUE, parallelAngle = 0.1,
               # vTrans = 210,
               vsize = vertex_size, 
               # esize = 7,
               asize = 12,
               mar = margin_size,
               label.cex = label_size,
               loopRotation = loop_rotation,
               minimum = 0, 
               maximum = maxedge,
               labels = c(# "Autonomy", "Competence", "Relatedness",
                 "Pleasure", "Interest", "Importance",
                 "Situation\nrequires", "Avoid guilt/\nanxiety", "Another\nwants"),
               label.scale = FALSE,
               color = viridis::viridis(6, begin = 0.5, direction = -1),
               # pie = pred_obj$tverrors[[tp]][, 3],
               directed = TRUE,
               title = "A)",
               title.cex = title_size) 

tp <- 25
qgraph::qgraph(t(tvvar_obj$wadj[, , 1, tp]), 
               layout = "circle",
               edge.color = t(tvvar_obj$edgecolor[, , 1, tp]),
               curveAll = FALSE,
               # vTrans = 210,
               vsize = vertex_size, 
               # esize = 7,
               asize = 12,
               mar = margin_size,
               label.cex = label_size,
               loopRotation = loop_rotation,
               minimum = 0, 
               maximum = maxedge,
               labels = c(# "Autonomy", "Competence", "Relatedness",
                 "Pleasure", "Interest", "Importance",
                 "Situation\nrequires", "Avoid guilt/\nanxiety", "Another\nwants"),
               label.scale = FALSE,
               color = viridis::viridis(6, begin = 0.5, direction = -1),
               # pie = pred_obj$tverrors[[tp]][, 3],
               directed = TRUE,
               title = "B)",
               title.cex = title_size) 

tp <- 45
qgraph::qgraph(t(tvvar_obj$wadj[, , 1, tp]), 
               layout = "circle",
               edge.color = t(tvvar_obj$edgecolor[, , 1, tp]),
               curveAll = FALSE,
               # vTrans = 210,
               vsize = vertex_size, 
               # esize = 7,
               asize = 12,
               mar = margin_size,
               label.cex = label_size,
               loopRotation = loop_rotation,
               minimum = 0, 
               maximum = maxedge,
               labels = c(# "Autonomy", "Competence", "Relatedness",
                 "Pleasure", "Interest", "Importance",
                 "Situation\nrequires", "Avoid guilt/\nanxiety", "Another\nwants"),
               label.scale = FALSE,
               color = viridis::viridis(6, begin = 0.5, direction = -1),
               # pie = pred_obj$tverrors[[tp]][, 3],
               directed = TRUE,
               title = "C)",
               title.cex = title_size) 

```

Idiographic science, which tries to unveil person-level processes, does not aim to go inductively from data to universal or statistical laws that hold in hypothetical infinitely large populations [@gaylesUtilityPersonspecificAnalyses2013; @piccirillo2019FoundationsIdiographicMethods]. Instead, it applies general principles, such as universal properties of complex systems, to study how specific individuals behave in their particular contexts. Answering more than half a century of calls to expand focus beyond outcomes to processes, new technology in data collection and analysis has now made the idiographic approach possible [@hamakerModelingBASDysregulation2016]. The basic solution is to not average individuals and then model the behaviour of the averages, but to first model individuals, and then aggregate those models to search for commonalities [@wrightPersonalizedModelsPsychopathology2020a]. Recent work has made use of methods such as ecological momentary assessment [e.g. @burkeEcologicalMomentaryAssessment2017] to gather intensive longitudinal data on behaviour and determinants from one or more individuals which can then be represented as time-series. In the case of smoking, analyses of such idiographic data have yielded individualized models that can predict behaviour with stunning accuracy [for some individuals at least; @fisherGeneratingAccuratePersonalized2019; @soysterInvolvingStakeholdersDesign2019].

Coming back to the notion of mechanisms; if the mechanisms happen within an individual, we need to study them at the appropriate level, that is, within-individual. However, when we study individual time series data, it becomes quickly obvious that the methods used in the conventional approach for studying group averages (e.g. pre-post measurements with a long time between them) leave us wanting. Figure \@ref(fig:sampling-rate-plot) illustrates that if insufficient within-individual time points are sampled, a deceptively linear picture of the process emerges (see also [@schiepekRealtimeMonitoringPsychotherapeutic2016], p. 3). The same logic applies if we are studying groups but cannot rely on the means being informative due to a lack of power (as demonstrated in [@carelloWhyNonlinearMethods2005]). 

```{r sampling-rate-plot, fig.width = 8.27, fig.height = 11.69 * (2/3), fig.cap = "One of the time series recorded by the participant featured in previous figure. Dots indicate answers to a visual analog scale question on their relatedness need satisfaction, as posited by self-determination theory (y-axis), measured on different time points (x-axis). A) Measuring three time points—representing conventional evaluation of baseline, post-intervention and a longer-term follow-up—shows a decreasing trend. B) Same measurement on slightly different days shows an opposite trend. C) Measuring 15 time points instead of three would have accommodated both observed \"trends\". D) New linear regression line (dashed) indicates stationarity. E) Including all of the 122 time points, a more complete picture of the dynamics emerges."}

emadata <- readr::read_csv("./data/EMA_data_Moti_P10.csv")
# emadata <- readr::read_tsv("../motivation-dynamics/data/data.tsv") %>%
#   dplyr::filter(User == "Moti_P10") %>% 
#   tidyr::spread(Field, Value)

emadata <- emadata %>% 
  dplyr::mutate(date = as.Date(dateTime),
                relatedness = relatedness / 50 * 100) %>% # scale to 1-100
  dplyr::mutate(timediff = c(NA, diff(dateTime))) %>% 
  dplyr::filter(timediff > 15) %>% 
  dplyr::select(-timediff) %>% 
  dplyr::mutate(time = row_number())


# Data frame containing only three points, to create regression lines in the plot:
data2 <- data.frame(y = c(emadata$relatedness[3], 
                          emadata$relatedness[55],
                          emadata$relatedness[122]),
                    x = c(3, 55, 122))

nonlinplot_a <- emadata %>% ggplot(aes(x = time, y = relatedness, group = 1)) +
  geom_point(aes(x = 2, y = emadata$relatedness[3]), color = "red", size = 2) +
  geom_point(aes(x = 54, y = emadata$relatedness[55]), color = "red", size = 2) +
  geom_point(aes(x = 121, y = emadata$relatedness[122]), color = "red", size = 2) +
  geom_smooth(data = data2, aes(x = x, y = y), 
              method = "lm", se = FALSE, color = "red", linetype = "dashed") +
  papaja::theme_apa() +
  scale_y_continuous(name = NULL) +
  coord_cartesian(xlim = c(0, nrow(emadata)), ylim = c(0, 100)) +
  labs(title = "A)") +
  xlab(NULL)

data3 <- data.frame(y = c(emadata$relatedness[2], 
                          emadata$relatedness[56],
                          emadata$relatedness[121]),
                    x = c(2, 56, 121))

nonlinplot_b <- emadata %>% ggplot(aes(x = time, y = relatedness, group = 1)) + # geom_point()
  geom_point(aes(x = 3, y = emadata$relatedness[3]), color = "red", size = 2) +
  geom_point(aes(x = 55, y = emadata$relatedness[55]), color = "red", size = 2) +
  geom_point(aes(x = 122, y = emadata$relatedness[122]), color = "red", size = 2) +
  geom_smooth(data = data2, aes(x = x, y = y), 
              method = "lm", se = FALSE, color = "red", linetype = "dashed") +
  geom_point(aes(x = 2, y = emadata$relatedness[2]), color = "blue", size = 2) +
  geom_point(aes(x = 56, y = emadata$relatedness[56]), color = "blue", size = 2) +
  geom_point(aes(x = 121, y = emadata$relatedness[121]), color = "blue", size = 2) +
  geom_smooth(data = data3, aes(x = x, y = y), 
              method = "lm", se = FALSE, color = "blue", linetype = "dashed") +
  papaja::theme_apa() +
  scale_y_continuous(name = NULL) +
  coord_cartesian(xlim = c(0, nrow(emadata)), ylim = c(0, 100)) +
  labs(title = "B)") +
  xlab(NULL)

nonlinplot_c <- emadata %>% ggplot(aes(x = time, y = relatedness, group = 1)) +
  # Plot every 9th point, as well as some selected ones in the end:
  geom_line(data = emadata[c(2, 3, 55, 56, 122, 121, 8, 18, 27, 31, 44, 67, 76, 90, 104), ], 
            aes(x = time, y = relatedness), color = "darkgrey", size = 0.75) +
  geom_point(data = emadata[c(2, 3, 55, 56, 122, 121, 8, 18, 27, 31, 44, 67, 76, 90, 104), ],
             aes(x = time, y = relatedness), color = "blue", size = 0.8) +
  geom_point(aes(x = 3, y = emadata$relatedness[3]), color = "red", size = 2) +
  geom_point(aes(x = 55, y = emadata$relatedness[55]), color = "red", size = 2) +
  geom_point(aes(x = 122, y = emadata$relatedness[122]), color = "red", size = 2) +
  geom_smooth(data = data2, aes(x = x, y = y), 
              method = "lm", se = FALSE, color = "red", linetype = "dashed") +
  geom_point(aes(x = 2, y = emadata$relatedness[2]), color = "blue", size = 2) +
  geom_point(aes(x = 56, y = emadata$relatedness[56]), color = "blue", size = 2) +
  geom_point(aes(x = 121, y = emadata$relatedness[121]), color = "blue", size = 2) +
  geom_smooth(data = data3, aes(x = x, y = y), 
              method = "lm", se = FALSE, color = "blue", linetype = "dashed") +
  papaja::theme_apa() + 
  scale_y_continuous(name = "Relatedness need satisfaction") +
  coord_cartesian(xlim = c(0, nrow(emadata)), ylim = c(0, 100)) +
  labs(title = "C)") +
  xlab(NULL)

nonlinplot_d <- emadata %>% ggplot(aes(x = time, y = relatedness, group = 1)) +
  geom_line(data = emadata[c(2, 3, 55, 56, 122, 121, 8, 18, 27, 31, 44, 67, 76, 90, 104), ], 
            aes(x = time, y = relatedness), color = "darkgrey", size = 0.75) +
  geom_point(data = emadata[c(2, 3, 55, 56, 122, 121, 8,  18, 27, 31, 44, 67, 76, 90, 104), ],
             aes(x = time, y = relatedness), color = "blue", size = 0.8) +
  geom_line(data = emadata[c(2, 3, 55, 56, 122, 121, 8, 18, 27, 31, 44, 67, 76, 90, 104), ],
              aes(x = time, y = relatedness),
              stat = "smooth", method = "lm", linetype = "dashed", color = "blue", alpha = 0.4) +
  papaja::theme_apa() + 
  scale_y_continuous(name = NULL) +
  coord_cartesian(xlim = c(0, nrow(emadata)), ylim = c(0, 100)) +
  # Red highlight circles:
  geom_point(data = emadata[c(2, 3, 55, 56, 122, 121), ],
             aes(x = time, y = relatedness), shape = 1, color = "red", size = 3.5) +
  labs(title = "D)") +
  xlab(NULL)

nonlinplot_e <- emadata %>% ggplot(aes(x = time, y = relatedness, group = 1)) +
  geom_line(color = "darkgrey", size = 0.25) +
  geom_point(color = "blue", size = 0.8) +
  geom_smooth(method = "lm", linetype = "blank") + # Produce confidence interval only
    geom_line(stat = "smooth", method = "lm", linetype = "dashed", color = "blue", alpha = 0.4) +
  scale_y_continuous(name = NULL) +
  scale_x_continuous(name = "Response number") +
  # Red highlight circles:
  geom_point(data = emadata[c(2, 3, 55, 56, 122, 121), ],
             aes(x = time, y = relatedness), shape = 1, color = "red", size = 3.5) +
  coord_cartesian(xlim = c(0, nrow(emadata)), ylim = c(0, 100)) + 
  papaja::theme_apa() + 
  labs(title = "E)")

# Patchwork syntax for composite plot
library(patchwork)

fullplot <- 
((nonlinplot_a | nonlinplot_b) / 
 (nonlinplot_c | nonlinplot_d) / 
    nonlinplot_e)

show(fullplot)
```
  
In sum, to study individual behaviour change, we need to not only collect intensive longitudinal data on the individual-level, but we must also consider the time evolution of the phenomenon and apply statistical analyses which can accurately model non-stationary data. In the health psychology context, Bolger and Zee [-@bolger2019HeterogeneityTemporalProcesses] argue, that not only temporal processes need to be considered, but also the heterogeneity therein. Consistent with the idiographic approach outlined above, every individual may exhibit idiosyncratic dynamics. As we will see next, the possibilities are vast outside of linear world. 

## Nonlinear dynamics

```{r covid}

n0 <- 10
r0 <- 0.20
r0_2 <- 0.19
t1 <- 30
t2 <- 45

cases1 <- n0*exp(r0*t1) %>% round(., digits = 0)
cases2 <- n0*exp(r0*t2) %>% round(., digits = 0)

cases2_reduction <- n0*exp(r0_2*t2) %>% round(., digits = 0)

reduction <- cases2 - cases2_reduction

```

As mentioned above, the linear methods traditionally used in psychology (e.g., multiple linear regression, ANOVA, and other cases of the general linear model) view psychological phenomena as following benign gradual changes over time. While sometimes useful as approximations, the assumptions of linear models are usually violated in practice [@siegenfeldIntroductionComplexSystems2019]. Furthermore, linear models may be invalid when ceiling or floor effects are present [@gonzalezNonlinearityComplexityLimited2010; @verboonApplyingGeneralizedLogistic2020], or under _hysterisis_, when the temporal direction of a relationship matters for its impact (e.g. prevention is important precisely because it takes more effort to exit the state of having a lifestyle disease, than to enter it) [@kelsoHakenKelsoBunzModel2008; @rooijModelingDynamicsRisky2013]. 

While a reliance on linear models simplifies the analytical approaches needed to explore relationships between variables, it does not contribute to our understanding of how the world works, as "most of everyday life is nonlinear" [@strogatzNonlinearDynamicsChaos2018, p. 9] and outside the physical sciences, nonlinear systems are "the rule, not the exception" [@maySimpleMathematicalModels1976, p. 467]. As an intuitive example, consider that falling from 10 meters is likely to kill you, but falling from one meter does not make you 1/10th dead – in fact, it makes you stronger [@talebAntifragileThingsThat2012; @talebAntifragilityMathematicalIdea2013]. Or that eating twice the size of a normal meal rarely results in twice the pleasure. 

Nonlinear dynamics, on the other hand, can be very useful, but unintuitive, to grasp, as the world discovered during the COVID-19 pandemic: An exponential growth starting from `r n0` cases with a growth rate of `r r0*100`% can lead to `r cases1` cases by day `r t1`, and `r cases2` cases by day `r t2` -- whereas a mere 1% reduction in the growth rate would have resulted in approximately `r ((reduction/1000) %>% round(., digits = 0)) * 1000` fewer cases by that time. Theories and methods to understand non-linear change phenomena in individuals can provide different types of answers than linear analyses. The most important factors in predicting behaviour change may not be the strength of a variable’s relationship with behaviour (e.g. regression weights), but rather the type of fluctuation that the variable exhibits in response to an intervention [e.g. so-called fractal, power-law, or 1/f noise; @almuradComplexityMatchingRestoring2018; @bakSelforganizedCriticalityExplanation1987; @delignieresFractalDynamicsSelfesteem2004], or how fast the dynamics recover after shocks [@vanordenLivingPinkIntentionality2011]. Another key insight is that, while we cannot usually predict what the value of the next observation will be, we can predict which system states are possible, and evaluate the risks and opportunities for intervention from there.

Polynomial regression is perhaps the most commonly used model when linearity is questioned. This method allows for identifying curves that may better fit data on the relationships between variables than a straight line [e.g. the convex---upward-curving---relationship demonstrated between intentions and behaviour between individuals; @chatzisarantisRelationshipPhysicalActivity2019], and can also be used to represent non-linear changes that occur over time. Polynomial regression models do not, however, adequately capture the essence of complex systems; nonlinear, irregular changes, periodic peaks and plateaus, and with recoveries after negative shocks and deterioration after positive ones [@hofmansKCentresFunctionalClustering2017, p. 2].

When we consider the situation where all components of a system interact, many features evident in everyday life but ambiguous in linear modelling become salient. Long periods with no discernible changes in outcomes might be followed by short bursts with large shifts. For example, a person’s conscious intention to smoke may remain stable, while social norms keep changing, until one day a seemingly innocuous event causes the person to quit. When a system finally reaches a "tipping point" (e.g. an individual’s behaviour changes), conventional analytic methods have difficulty determining whether the effect was caused by a critically important incident, or by less obvious, small, cumulative effects over time which preceded the so-called _phase transition_. Obviously, in such situations, the consequences of an incident (i.e.  the camel’s back breaking) do not relate linearly to the intensity of the event (i.e. loading the last straw on the camel). This is a common dynamic in complex systems [@talebBlackSwanCairo2011], but it is extremely difficult to evaluate if information regarding the system is only available for a few points in time. Intensive longitudinal data is therefore needed.  

# Empirical solutions

To model intensive longitudinal data, models developed within the literature on time series analysis are necessary [@bradleyNonlinearTimeseriesAnalysis2015; @wrightPersonalizedModelsPsychopathology2020a]. A time series in this case is a sequence of values representing one variable in one individual, and time series analysis consists of methods for studying time evolution of one or more data generating processes. 

The most common modelling framework, lag-1 autoregression, uses one previous time point as input to predicting the next one. In behavioural science, vector autoregression---vectors being sequences of numbers, representing values of variables---is often used to test the effects of several variables on the outcome of interest. One drawback of such autoregressive models is that they assume that there exists an average value around which the process fluctuates, which also motivates the common practice of "detrending". In detrending, the researcher transforms the data by fitting a linear regression line and continuing the analysis with the residuals, often not taking into account that there can be several trends in subsections of the data (i.e. the trend is non-stationary), which all contribute to what the linear model interprets as normally distributed "errors". Moreover, the supposed mean value---as well as variance around it---may not remain the same across time (i.e. the level is stationary), and the impact of previous time points on future ones is assumed to remain constant [@bringmannChangingDynamicsTimevarying2017, p. 5]. One way to overcome this particular shortcoming, is to let the parameters in autoregressive models vary across time, leading to the time-varying autoregressive model depicted in Figure \@ref(fig:tv-var). But even time-varying autoregressive models operate under the linear regression framework, with its accompanying assumptions, such as normally distributed errors. Furthermore, in Figure \@ref(fig:tv-var) we have limited ourselves to investigating the lag-1 relationships, whereas long-range dependencies are common in ecological momentary assessment data [@delignieresFractalDynamicsSelfesteem2004; @hasselmanStudyingComplexAdaptive2020; @olthofPsychologicalDynamicsAre2020; @olthof2020ComplexityPsychologicalSelfratings]. 

Regression-based approaches including time-varying autoregressive models are only appropriate when the dynamics of all variables in the model conform to the required assumptions. Empirical researchers have a wide variety of assumption tests at their disposal. The supplementary website (section https://git.io/JfLmm) presents a plethora of these tests applied to a sample of 20 individuals collecting motivation data for nine variables. We can see that many or most time series indeed seem to exhibit non-stationary trends and levels, as well as non-linearities. Also, longer time series reject more of the assumptions, as the deviations from assumptions are not necessarily present in small samples, and larger samples confer higher statistical power. This does not suggest that we ought to only gather short time series, as doing so would limit our abilities to detect deviations from assumptions and generalise to data outside the sample.

There are many ways to study nonlinear change processes in complex systems. Behavioural researchers may find the generalised logistic model [@verboonApplyingGeneralizedLogistic2020] a good starting point. This method produces readily-interpretable parameters indicating the floors and ceilings of the variables intervened upon, as well as the growth rate and timing of changes. Researchers may also be interested in identifying critical transformations taking place in a system (e.g. a person’s motivational system). In complex systems, these shifts may be preceded by warning signs such as increased turbulence (quantified as e.g. _dynamic complexity_; [@schiepekIdentificationCriticalFluctuations2010]), or critical slowing down (i.e. heightened autocorrelations in a time series), before (re)lapses occur [@leemputCriticalSlowingEarly2014; @wichersCriticalSlowingPersonalized2016]. In clinical psychology interventions, intensive monitoring of psychopathological symptoms has allowed researchers to examine symptoms’ variability, autocorrelations and other indicators of dynamics. This has yielded considerable advances in the prediction of phase transitions between adaptive and maladaptive states during interventions in the field of psychopathology [@jeronimusDynamicSystemPerspectives2019; @olthofCriticalFluctuationsEarlyWarning2019; @olthofDestabilizationSelfratingsPsychotherapeutic2019; @schiepek2020ConvergentValidationMethods]. A conceptual replication was recently done in a population undergoing a weight loss intervention, where [@chevance2020IdiographicDaytodayFluctuations] found that sudden drops in physical activity levels could be predicted by the emergence of erratic fluctuations in day-to-day activity. While @olthofDestabilizationSelfratingsPsychotherapeutic2019 identified that the presence of critical fluctuations was a key indicator of the effectiveness of psychotherapy for mood disorders, this has not yet been investigated in other behaviour change contexts. 

In the next section, we exemplify one particular family of analysis methods, recurrence quantification, due to its suitability for analysing many existing longitudinal data sets while making fewer a-priori assumptions. This enables us to observe more granularity in the dynamics, than allowed by e.g. multilevel models, which treat individuals as departing from group-level means according to a known distribution [@piccirillo2019FoundationsIdiographicMethods].  

## Modeling complex time series data with Recurrence-based Analyses

To explore the dynamics of a phenomenon while making no assumptions about distributional shapes of observations or their errors, about linearity, or about the time-lags involved, researchers can perform Recurrence Quantification Analysis, which provides a robust visual intuition about the organisation of a system (recall from Table 1 that in complex systems, the organisation of components can be more important than the components themselves). An in-depth walkthrough of the analysis with code is provided in the supplementary website (see section https://git.io/JfLs3), hence we will be brief in the background and focus on the results. The data we use for this demonstration is a single participant, who at each time point was prompted by a mobile app to answer six questions about their motivation. A more detailed exposition of the data is found at the supplementary website (see section https://git.io/JfLmQ). 

Recurrence networks display relationships between multivariate observations in a time series in an intuitive way, which in the case of multidimensional Recurrence Quantification Analysis (see [@wallotAnalyzingMultivariateDynamics2018] for a tutorial) can be thought of as displaying a type of multivariate "correlation", indicating which occasions (be they single values or combinations of values of different variables as a system state) repeat a particular previously-observed pattern. These patterns or configurations can be thought of as being produced around attractors, towards which the system is drawn. 

Figure \@ref(fig:recnet) demonstrates such a multidimensional recurrence network, where each point is a measurement occasion, and recurrences of "motivation profiles" (consisting of the six motivation-related questions) are shown as links between time points. We can see that most of the recurrences take place in the second half of the data. In addition, all the patterns (or, equivalently: attractors, profiles, configurations) which occur only once, take place in the first half of data collection. Having only measured the first 50% of observations would have missed much of the action.

```{r recnet, fig.width = 8.27/2, fig.height = 8.27/2, fig.cap = "Weighted multidimensional recurrence network. Each circle (\"node\") is a measurement occasion, numbers indicate their running number, and colors represent different motivation profiles. These profiles are configurations of six variables, and can be conceived of as attractors. Lines indicate the same motivational state reoccurring at a later time point. Yellow nodes indicate configurations connecting to that with the highest strength centrality (i.e. number of connections weighted by the similarity of the connected nodes), red nodes connect to the second strongest which is not connected to the strongest, followed by purple and blue. Grey nodes depict uncategorised configurations which occur at least twice, and white ones the configurations, which only occur once. Nodes that are larger, are connected to more other nodes. Drawn with R package casnet [@hasselman2020CasnetToolboxStudying]."}

# knitr::include_graphics(path = "./figures/rqa_multiplot.png")

grid::grid.raster(png::readPNG("./figures/recnetwork.png"))

```

Figure \@ref(fig:attractors) shows the participant's four major repeating profiles, which correspond to the coloured nodes in Figure \@ref(fig:recnet). We can observe that about a fifth of the participant's responses fall into a relatively balanced profile, while ca. 15% indicates what self-determination theory [@ryanSelfdeterminationTheoryBasic2017] would consider an "optimal" motivation profile -- high in autonomous forms of motivation and low in controlled ones.

```{r attractors, fig.width = 11.69/2, fig.cap = "Main profiles corresponding to the colors indicated in the previous plot. See supplementary website (section https://git.io/JfLmS) for a thorough exposition."}

# knitr::include_graphics(path = "./figures/rqa_multiplot.png")

grid::grid.raster(png::readPNG("./figures/attractors.png"))

```

If we want to simplify, the information in Figure \@ref(fig:recnet) can be collapsed to the relative frequencies of a particular state following another, that is, demonstrating a _transition network_. It answers the question "If you are in state x, what is the probability of transitioning to state y?". Figure \@ref(fig:transition) (panel A) depicts the relative frequencies on which a particular state precedes the next -- in the case of linear dynamics, the matrix would be symmetric with similar values observed above and below the diagonal drawn from bottom left to the top right cell. Panel B of Figure \@ref(fig:transition) presents the same information as a network.

```{r transition, fig.width = 11.69/2 + 2, fig.cap = "A transition network.  If the system is in the configuration labelled 1st, based on the relative frequencies of observed transitions, there is 38% chance it stays in the same configuration, and a 16% chance it transitions to one of the uncategorised states – that is, states that are less strong than the state labelled 4th, but which appear more than once. Note that the .00 values in panel A do not signify this transition is impossible, only that it did not appear once during the data collection period."}

grid::grid.raster(png::readPNG("./figures/transition_grid_network.png"))

```

To distinguish whether the results reflect non-linear structure in the data or are merely a product of randomness, the researcher can take advantage of a technique called surrogate data analysis [@schreiber2000SurrogateTimeSeries]. The analysis is presented in the supplementary website (section https://git.io/JqRTQ), but in brief, temporally disordered versions of the data---called “surrogates”---are created, and the observed data is compared to those. The surrogates represent the hypothesis that the data were generated by a rescaled Gaussian linear process. This means that, by analysing the surrogates, we ask whether the data can be understood to have arisen from a process, that is essentially stochastic and linear instead of highly interdependent and non-linear. The analysis indicates, that it would indeed be very unlikely to see these results, if the dynamics were Gaussian..

Following the idiographic approach outlined earlier, this information could now be used to develop a personalised intervention. Designer of such an intervention -- possibly the person themself, using self-enactable behaviour change techniques [@knittleCompendiumSelfenactableTechniques2020] -- has several new perspectives to consider: How could transitions to the "optimal" profiles (2nd and 4th) be increased? Could they be made more "sticky", and the 1st state less so (as indicated by self-loops in Figure \@ref(fig:transition), panel B)? Could e.g. just-in-time interventions -- such as prompts on a mobile device -- be used to make the person more aware of the state they are in, and remind them of techniques that are helpful to staying in, or leaving, that state?

To summarise: Having looked at all the time scales instead of just the previous time points, while not restricting ourselves to linear dynamics, we observe the features outlined in Table 1. The repeating states are connected across time, demonstrating interconnectedness instead of independence. They are not equally spaced in time, demonstrating non-stationarity and hence non-ergodicity. In addition, going from state a to state b does not generally happen with the same probability as the reverse, hence demonstrating non-linear dynamics. These features would by default be smoothed out in regular approaches, but they are potent information for formulating interventions and understanding the dynamics at play.

# Discussion 

Applied behavioural sciences have always studied phenomena, like behaviour change mechanisms, which take place within complex ecological systems [@bronfenbrennerEcologyDevelopmentalProcesses1998], but in the majority of cases we have tried to understand these phenomena using linear models, when the tools of complexity science would have been more appropriate [@navarroTakingTimeSeriously2015]. Behavioural scientists have an opportune moment to start considering complexity, as the field of behavioural intervention research is now taking committed first steps in this direction [@craigTakingAccountContext2018; @skivingtonDevelopingEvaluatingComplex2018], and there is a growing interest toward intervention programme theories that explicitly model complex aspects, such as recursive causality, disproportionate relationships, "tipping points", and emergent outcomes [e.g. @rogersUsingProgrammeTheory2008a]. In addition, analytical methods that are compatible with complexity science, have recently been, and are increasingly being, developed [@hasselman2020CasnetToolboxStudying]. 

In an empirical dataset of 20 individuals, we have shown similar results to a plethora of studies, which find complex dynamics in ecological momentary assessment data [e.g. @delignieresFractalDynamicsSelfesteem2004; @navarroApproachStudyDynamics2007; @olthof2020ComplexityPsychologicalSelfratings; @olthofPsychologicalDynamicsAre2020]: Nonlinear, non-ergodic, non-independent dynamics, which defy traditional assumptions. Moreover, the empirical case example of a single person shows the importance of observing change over long periods of time: There are several reoccurring "ways of being", patterns of experience which change dynamically within the individual. These patterns of recurring states, and the nature of which states recur, are likely to be different across individuals. An added benefit of using multivariate recurrence-based methods [@hasselmanStudyingComplexAdaptive2020] is that they do not require partialing out variance, allowing the researcher to examine conceptually overlapping variables from multiple theories in the same analysis [@gordon1968IssuesMultipleRegression]. This makes intervention process evaluation possible from a more holistic perspective -- that is, looking for changes in e.g. attractors or complexity measures.

Critically appraising the often hidden assumptions of models, especially in the context of complex systems such as human behaviour change interventions, is necessary for understanding the phenomena of interest and building a credible science. While researchers who study stable phenomena and only wish to draw group-level inferences (e.g. to select promising public health interventions) are probably best served with traditional models, this is rarely the case for psychologists and behaviour change intervention researchers who wish to understand how behaviour changes. For theory to advance, assumptions need to be justified: We cannot conclude both that our models for empirical testing omit crucial facets of reality, and at the same time imply real-life consequences. We propose that a more fruitful approach would be to model coupled processes with individual-level psychological data from intensive longitudinal designs, using analyses which are reasonably free from assumptions regarding independence, ergodicity and linearity. By studying what other sciences know about change processes in complex systems, and replicating studies where the ideas have been applied to human behaviour change, researchers can work towards uncovering more general principles of behaviour change. As Molenaar [@molenaarPsychologicalMethodologyWill2007, pp. 38-39] pointed out, "the set of person-specific time series models thus obtained then can in the next step be subjected to standard analysis of inter-individual variation in order to detect subsets of subjects who are homogeneous with respect to particular aspects of the dynamical laws concerned". In other words, information obtained from individual-level studies of dynamic patterns can then possibly inform models of larger groups, leading to better (or at least humbler and more nuanced) social scientific theories [@smaldinoTheoryDevelopmentAgentBased2014]. Generating theory in this way would answer calls to address the issue of time more clearly in theories of health behaviour [@scholz2019ItTimeThink]. It could also lay the foundation for more formal theories of behaviour change to be developed [@chevance2020InnovativeMethodsPredicting],  as these typically hypothesise how relationships between variables unfold over time, and a more coherent correspondence between theoretical cycles and empirical cycles in behaviour change research [@vanrooij2021TheoryTestHow]. 

# Limitations

The field of complexity science and aligned novel methods is fast-moving, with new developments always on the horizon. However, there remain many practical and methodological barriers to fully embracing the complexity perspective in behaviour change research. Many of these barriers relate to data collection. While the development of smartphones and an array of other devices for ambulatory assessment allow the convenient collection of intensive longitudinal data, there are few stable and user-friendly open source options. This has resulted in large variability in the data collection tools used to produce intensive longitudinal data [@trullAmbulatoryAssessmentPsychopathology2020]. Ensuring good adherence to these forms of data collection can be a challenge for researchers. For participants, adapting to intensive assessment is a behaviour change in itself -- particularly if they are required to use a specific device or smartphone application. Although measurement burst designs [@stawski2015MeasurementBurstDesign] might potentially alleviate some challenges, they bring about other, perhaps graver ones, such as a mismatch between the sampling rate and the time scale on which the phenomenon in interest unfolds.

Long time series can be time-consuming and effortful to collect. It also creates a much greater burden on participants than traditional questionnaires and few timepoints only. However, in behaviour change research and health psychology, much of the core research interests of our theories---influences on behaviours---have traditionally been subjective factors (e.g., sense of self-efficacy, motivations and motives, outcome expectancies), only---by definition---accessible via self-report. This presents an undeniable practical challenge, along with the fact that observations generally need to be spread equidistantly in time, allowing one individual to collect only 1-2 data points per day. Still, examples of more than a hundred time points being collected are found from weight loss maintenance [e.g. @kwasnickaNof1StudyWeight2017] to psychotherapy [e.g. @olthofCriticalFluctuationsEarlyWarning2019], with some studies collecting more than a 1000 observations [@delignieresFractalDynamicsSelfesteem2004]. 

A number of methodological challenges for the study of dynamic systems in behavioural science have been identified [@hamakerNoTimePresent2017], including measurement reactivity, the optimal choice of measurement intervals, and measurement quality. To properly address measurement reactivity, it is necessary to know whether the anticipation of measurement or the self-monitoring process itself (or both) interact with the outcomes of interest. Choosing an optimal measurement interval requires knowing the timescale of the behaviour change dynamics, which is rarely well understood. As regards measurement quality, we still lack a comprehensive approach to developing and establishing the quality of momentary measures of psychological constructs. Ensuring the validity and reliability of these measures can be difficult due to the requirement to use few items, not to mention that the questionnaire scales are themselves bounded, whereas experience hardly is. One solution for this is to inspect change profiles of responses [@hasselmanStudyingComplexAdaptive2020] instead of raw scores. Another solution would naturally be tapping into wearable data; for example, electronically activated recorders [@kaplan2020BestPracticesElectronically] are maturing as a technology, and complexity methods have already been applied to  physical activity data during a weight loss intervention [@chevance2020IdiographicDaytodayFluctuations]. <!-- This alleviates, to an extent, the possibility of sudden shifts of how a person conceives the relation between their internal states and the slider or ordinal scale to which they ought to be projected to, as well as the likely situation where each answer is made relative to previous answers. Following the analogy of Taleb, Canetti, Kinda, Loukoianova & Schmieder -@talebNewHeuristicMeasure2012; "Using an inaccurate tape measure will give a false reading of a child’s height [...] However, if one uses the same tape measure over time, it will give a reliable test of whether the child is growing". -->

# Conclusion

When a study finds that variables have explained an unsatisfactory proportion of behaviour, researchers often follow the pattern seen in social and organisational sciences and conclude that either: "(a) significant, explanatory variables have been omitted from the study, (b) the measurement instrument is too imprecise and 'rough', or that (c) the random or stochastic part of the problem has overwhelmed the patterned part" [@mathewsWhyStudyComplexity1999, p. 453]. But if the result stems from a statistical model that makes unfounded assumptions regarding independence, ergodicity and linearity, is it any wonder that it fails to satisfactorily describe reality? In this paper, we have attempted to show that many common modelling strategies fail to adequately capture real-world dynamics of behaviour change, and that a change in approach can advance our understanding of behaviour and behaviour change processes. Behaviour change researchers should further utilize intensive longitudinal designs to collect individual-level psychological and behavioural data from participants, and should increasingly analyse such data with methods that are reasonably free from assumptions of independence, ergodicity and linearity. This has practical implications from replicability to outcome and intervention selection. In our view, further embracing complexity science and its methods will advance research on behaviour change and could unearth new evidence of the dynamics of behavioural processes.

\noindent 

## Declarations  

### Conflict of interest

The authors report no conflict of interest.

### Ethical approval

The research proposal was reviewed by the University of Helsinki Ethical Review Board in the Humanities and Social Sciences.

### Authors' contributions

MH wrote the analysis code, including the full online supplement, formulated the initial draft of the manuscript and revised it. KK, CN, NH and FH revised the paper for its intellectual content in several rounds. MH, KK and NH contributed to collecting the empirical data. FH wrote and adapted the R package used for analysis, and provided expertise in the analytical approach. All authors read and approved the final manuscript.

### Funding

MH and KK were supported by Academy of Finland (grant number 295765 to principal investigator NH). NH was supported by an Academy of Finland Research Fellowship (grant number 285283). 

### Reporting 

We report all data exclusions, all manipulations, and all measures in the study on the supplementary website (section https://git.io/JfLmQ). 

### Data, materials, and online resources

The Open Science Framework repository for this project is at https://osf.io/hnpse/. This manuscript was created with the R package papaja [@aust2020PapajaPreparingAPA], and the source code is available at https://git.io/JfLY1. All analyses and code are available on the supplementary website (https://git.io/Jffrm, GitHub repository at https://git.io/JfLYr). Data is available at https://git.io/JfLYo.

### Acknowledgements

We would like to thank Daniele Proverbio and Jonas Haslbeck for their comments, as well as Anton Kunnari for the idea of transition networks.

\newpage

# References

\begingroup
\setlength{\parindent}{-0.5in}
\setlength{\leftskip}{0.5in}

<div id = "refs"></div>
\endgroup