Clinical factors and early life experiences associated with therapeutic alliance development in treatment for depression or binge eating

Abstract

Objective:

This study examines childhood and clinical factors theorized to impact therapeutic alliance development over the course of psychotherapy.

Method:

Raters assessed the therapeutic alliance of 212 client-therapist dyads, participating in two randomized controlled trials of schema therapy and cognitive behavioural therapy for binge eating or major depression, at three time points. Linear mixed models were used to characterize therapeutic alliance development over time and assess the influence of childhood trauma, perceived parental bonding, diagnosis and therapy type on scores.

Results:

Participants differed in initial alliance ratings for all subscales but had similar growth trajectories in all but the patient hostility subscale. A diagnosis of bulimia nervosa or binge eating disorder predicted greater initial levels of client distress, client dependency and overall client contribution to a strong therapeutic alliance, compared with a diagnosis of depression. Therapy type, childhood trauma and perceived parental bonds did not predict alliance scores.

Conclusion:

Findings highlight the potential influence of clinical and personal characteristics on alliance strength and development, with implications for maximizing treatment outcomes through anticipating and responding to these challenges.

Keywords:

• therapeutic alliance
• linear mixed models
• childhood trauma
• parental bonding

Clinical or methodological significance of this article: A moderate association exists between therapeutic alliance quality and psychotherapy outcomes, however little is known about specific client characteristics that predict alliance strength. This article identifies specific diagnoses as predictors of client contribution to the therapeutic alliance, as well as early childhood variables of interest for future alliance research. The article also illustrates potential problems that may arise with therapeutic alliance measurement in research, such as reduced variation in scores due to the use of independent raters or controlled therapy environments.

The therapeutic alliance is commonly defined as the strength of emotional bond shared by the client and their clinician, in addition to agreement on long-term therapy goals and the tasks undertaken to achieve them (Ardito & Rabellino, 2011; Bordin, 1979; Flückiger et al., 2018). It is now well established that a strong therapeutic alliance is linked with favourable outcomes in any form of treatment (Flückiger et al., 2018) and that certain therapist attributes and skills are in turn linked with a better quality of therapeutic alliance (e.g., Ackerman & Hilsenroth, 2003; Nienhuis et al., 2018). How client characteristics contribute to the therapeutic alliance is less clear. Identifying client factors that may foster or diminish therapeutic alliance strength is an important step in creating targeted approaches to alliance building.

The broad interpersonal patterns of clients have proven to be more consistent predictors of therapeutic alliance, compared to clinical or demographic factors (Cheng & Lo, 2018; Constantino et al., 2005; Garner et al., 2008; Saunders, 2001). For example, attachment style (Bowlby, 1988) has been associated with alliance strength in several studies (Diener & Monroe, 2011; Folke et al., 2016; Smith et al., 2010). This makes intuitive sense, as the manner an individual habitually perceives and responds to others in everyday life is likely to carry into the therapeutic context and directly influence the relationship formed with a therapist. Additional early experiences that may strongly influence later interpersonal styles, such as perceived bonds with caregivers and adverse early experiences, are therefore of interest to therapeutic alliance research. Client ratings of childhood bonds with parents have been found to be associated with therapeutic alliance quality in a number of studies (Bovard-Johns et al., 2015; Hersoug et al., 2002, 2009), however it is unclear whether specific aspects of these relationships (such as care or protection) are more closely associated with therapeutic alliance in later therapy (Mallinckrodt, 1991; Mallinckrodt et al., 1995) or whether there are differing associations with paternal versus maternal bonds (Bovard-Johns et al., 2015; Mallinckrodt et al., 1995). Similarly, an association between adverse childhood events and later therapeutic alliance quality has been identified, albeit inconsistently (Bovard-Johns et al., 2015; Eltz et al., 1995; Keller et al., 2010; Lubman et al., 2007; Lysaker et al., 2011; Paivio & Patterson, 1999; Reynolds et al., 2017). Further exploration of potential associations between these early experiences and therapeutic alliance is needed in order to better understand and respond to potential barriers in therapeutic alliance development, and thus to therapy success.

The current study also contributes to the existing body of research on client factors that may influence therapeutic alliance quality by addressing common methodological issues, which may explain the inconsistencies in the research noted above. The broad, subjective and dynamic nature of the therapeutic alliance construct itself has led to queries about its validity in recent years. For example, in a recent review Horvath (2018) identified over 70 instruments with varying degrees of scientific validation in use as measures of therapeutic alliance. There is also the question of from whose perspective the therapeutic alliance is best measured: the client’s, the therapist’s, or an independent observer’s? Correlations between client and therapist ratings of alliance quality have typically been low, with one meta-analysis finding a mean correlation of only r = .36 between the two rater types (Tryon et al., 2007). Client ratings are utilized more often in research as they are typically found to be the most predictive of treatment outcomes (Tschuschke et al., 2020), however social desirability and halo effects remain a concern (Cecero et al., 2001; Elvins & Green, 2008). Independent raters have thus been proposed as a more reliable means of capturing therapeutic alliance strength (Fenton et al., 2001).

Further, the quality of any therapeutic alliance is constantly in flux, both during and between sessions, yet is often measured at one time point or represented as an averaged score. The use of more complex research designs and analyses able to capture the dynamic nature of the alliance has resulted in a newly proposed model of therapeutic alliance, consisting of two components that differentially influence the therapeutic change process (Zilcha-Mano, 2017; Zilcha-Mano & Fisher, 2022). The first is a trait-like form of alliance. Clients with sufficient underlying capacity to form satisfying and productive relationships in general are expected to form a strong alliance with a “good enough” therapist early on in treatment, which then fosters engagement with the technical interventions of therapy (e.g., cognitive restructuring). A growing number of studies suggests these clients have the best therapy outcomes overall (Zilcha-Mano & Fisher, 2022). In contrast, clients with a more fraught relationship history may form unsteady alliances with their therapists, at least initially, and thus struggle with the suggested therapy interventions. The second component of the alliance is state-like, referring to changes in alliance strength across therapy. Significant positive changes in state-like alliance, such as those occurring after a particularly powerful use of empathy or successful rupture repair, appear to serve a directly curative function for clients with low trait-like alliance in particular (Zilcha-Mano & Fisher, 2022). For example, in a recent study of 42 participants undergoing emotion-focused therapy, changes in alliance strength were more closely linked with treatment outcome for participants who initially struggled to form a strong alliance (Harrington et al., 2021). Distinguishing between trait- and state-like alliance, for example by disentangling between and within subjects changes using linear mixed modelling (Lorenzo-Luaces & DeRubeis, 2018), is now considered important for reconciling previously contradictory research findings (Zilcha-Mano & Fisher, 2022) and developing a greater understanding of how context contributes to alliance strength and outcomes, with important practical implications for individual therapy.

Objective

To address these common methodological issues in alliance research, therapeutic alliance in the current study was measured by independent raters at three time-points throughout treatment, using two well-validated instruments completed from an observer perspective. Linear mixed models were used to assess between-subjects differences in alliance scores at baseline (to approximate trait-like alliance), averaged alliance changes over time, and within-subjects alliance changes over time (to approximate state-like alliance changes), and to examine whether factors such as therapy type, diagnosis, childhood trauma and perceived parental bonding before the age of 16 accounted for any observed variation. It was hypothesized that a poorer initial quality of alliance would be observed for participants who reported having experienced one or multiple traumatic childhood events in childhood, or who perceived their parental bonds to have been low in care, high in overprotection or both. A dose–response association, wherein more pathological parenting behaviours or traumas of an interpersonal nature would be associated with weaker initial therapeutic alliance, was also expected. The sample for the current study was drawn from two completed randomized controlled trails (RCTs) that compared schema therapy and cognitive behavioural therapy (CBT), and it was hypothesized that state-like changes in therapeutic alliance would also be more pronounced in schema therapy versus CBT, due to the explicit alliance focus of schema therapy. No specific hypotheses were made about possible differential effects of participant primary diagnosis, or parental bonds with mothers versus fathers or the bond dimensions of care versus overprotection on therapeutic alliance formation. Their influence on the therapeutic alliance was examined in an exploratory manner.

Methods

Sample

Data for the current study were derived from two RCTs (N = 212) conducted in Christchurch, New Zealand. The first compared schema therapy (Young et al., 2003) and CBT (Beck, 2011) for individuals diagnosed with major depression (N = 100, mean age 38 years; Carter et al., 2013), henceforth referred to as the “major depression (MD)” study. The second compared schema therapy, standard CBT (Fairburn et al., 1993), and appetite-focused CBT (McIntosh et al., 2007) in the treatment of individuals diagnosed with bulimia nervosa or binge eating disorder (N = 112, mean age 35 years; McIntosh et al., 2016). This study is henceforth referred to as the binge eating (BE) study, as binge eating is a shared symptom of both diagnoses. Participants were mostly female (N = 181) and self-identified as New Zealand European (MD: 85%, BE: 67%), followed by other Caucasian (MD: 8%, BE: 17%), Māori (MD: 4%, BE: 10%), Asian (MD: 1%, BE: 4%) and other non-Caucasian ethnicities (MD: 2%, BE: 2%). Therapy was delivered by six experienced female clinical psychologists in the MD study, and four female clinical psychologists in the BE study. Further details of the two samples, such as exclusion criteria, clinical characteristics, participant flow data and results comparing treatment outcomes and client predictors of response have been published elsewhere (Carter et al., 2013, 2018; Jordan et al., 2014; McIntosh et al., 2016).

Measures

Independent variables

Primary psychiatric diagnosis and comorbidity were assessed by trained clinicians at baseline using the research version, patient edition of the Structured Clinical Interview for DSM-IV-TR (SCID-I/P; First et al., 2002). Perceived parental bonding before the age of 16 was measured by the Parental Bonding Instrument (PBI; Parker et al., 1979). The PBI lists 25 behaviours and attitudes indicating care (warmth/attention vs. coldness/rejection) and overprotection (autonomy/independence vs. over-control/intrusiveness), on which each caregiver is retrospectively rated using a Likert-type scale from 0 (very unlike) to 3 (very like). Subscales can be dichotomized into “high” and “low” categories, with high care defined as scores above 27 for female caregivers and above 24 for male caregivers and high overprotection defined as scores above 13.5 for female caregivers and above 12.5 for male caregivers. The PBI has demonstrated high construct and concurrent validity and long-term test-retest reliability (Parker, 1983; Wilhelm et al., 2005). Traumatic events experienced in childhood were assessed as part of a larger clinician interview in which participants were asked about witnessing or experiencing disasters, accidents or war (i.e., impersonal trauma) or sexual, psychological or physical forms of abuse before the age of 16. Questions were based on those used to assess retrospective experiences of physical and sexual abuse in the longitudinal Christchurch Health and Development Study (Fergusson et al., 1996).

Dependent variables

Ratings of therapeutic alliance quality were made of audio-recorded therapy sessions, using two observer measures: The Vanderbilt Psychotherapy Process Scale (VPPS; O'Malley et al., 1983) and the Vanderbilt Therapeutic Alliance Scale (VTAS; Hartley & Strupp, 1983). Both assess aspects of the therapy relationship that have been theoretically linked with Bordin’s three dimensions of therapeutic alliance. The VPPS contains 80 items assessing the demeanour and behaviour of therapists and clients in a full therapy session, measured on a 5-point Likert-type scale from 1 (not at all) to 5 (a great deal). These items form eight subscales, labelled patient participation, patient hostility, patient exploration, patient distress and patient dependency, therapist exploration, therapist negative attitude, and therapist warmth (the term client will be used synonymously with the term patient in the current study, to retain the original subscale names and reflect preferred terminology). The modified version of the VTAS contains 37 items, which sum to create three subscales (Krupnick et al., 1996). Two of these subscales correspond to factors identified via factor analysis (Krupnick et al., 1996) and were used in the current study. The therapist factor subscale describes therapist behaviours that facilitate a strong therapeutic alliance (e.g., warmth, empathy), while the patient factor subscale measures alliance-promoting behaviours of the client (e.g., making an effort to carry out therapeutic interventions, openness, lack of defensiveness) and of their interactions with therapists (e.g., showing enthusiasm together, sharing a common viewpoint, acceptance of roles and responsibilities). All items measure the extent that these behaviours and attitudes were present, scored on a 6-point Likert-type scale ranging from 0 (not at all) to 5 (a great deal). Statistical indicators of good inter-rater reliability, internal consistency and predictive validity have been demonstrated for the VPPS (Suh et al., 1989; Windholz & Silberschatz, 1988) and the VTAS (Hartley & Strupp, 1983; Krupnick et al., 1996).

Procedure

Participants were recruited via advertisement, self-referral or referral from health professionals. Baseline assessments included self-report questionnaires (including the PBI), structured clinical interviews (including questions about childhood trauma), physical and neuropsychological assessment. Participants were then randomized equally between interventions in the MD study, while the BE study had a three-arm parallel group randomization. Therapies were delivered weekly for six months, followed by six months of monthly sessions. Therapeutic alliance ratings were made by trained observers (five graduate students in the BE study and seven in the MD study), who were naïve to treatment condition. One rating per participant was completed for early (sessions 2–5), middle (with exact session numbers differing among participants due to flexibility in the total number of sessions) and late (penultimate four sessions) stages of treatment, except for clients who ended treatment prematurely. At least twenty percent of rated sessions were randomly selected for rating by a second rater to evaluate inter-rater reliability.

Statistical Analysis

Inter-rater reliability was assessed by the one-way random intraclass correlation coefficient (ICC; Shrout & Fleiss, 1979), with estimates between 0.5 and 0.7 considered to indicate moderate reliability, above 0.7 good reliability and above 0.9 excellent reliability (Koo & Li, 2016). Inter-rater agreement was assessed using coefficients of variation (CV; a unit-free measure of dispersion representing the ratio of the standard deviation to the mean) with scores below 0.25 indicating high agreement. ICCs and CVs were calculated in SPSS (IBM Corporation, 2020), using 67 pairs of ratings from the BE study and 64 pairs from the MD study. Internal consistency of primary ratings made in the first phase of therapy was assessed using coefficient alpha (Cronbach, 1951) and coefficient omega (McDonald, 1999), with estimates of 0.7–0.9 accepted as good (Tavakol & Dennick, 2011). The traditional use of coefficient alpha to assess internal consistency has come under increasing criticism due to the restrictive assumptions of the statistic, thus coefficient omega is included for comparison purposes (McNeish, 2018). ICCs, CVs and alpha coefficients were calculated using SPSS (IBM Corporation, 2020), and omega coefficients using the “Psych” package in R (Revelle, 2019).

Linear mixed models were fitted using the “nlme” package in R (Pinheiro & Bates, 2022). Mixed models allow accurate analysis of hierarchical data (e.g., participants nested within therapy groups or observations nested within participants). A null model was first fitted to calculate ICCs (LeBreton & Senter, 2008; Shrout & Fleiss, 1979), with higher values representing greater variance in outcome attributable to the clustered nature of data. Several possible unconditional (i.e., without predictor) models were then compared in a stepwise fashion, using a combination of model fit indices, practical considerations and theoretical knowledge bases. Smaller Akaike Information Criterion (AIC; a conservative version of the −2 Log Likelihood that penalizes for model complexity) values indicated better fit when comparing models regardless of the absolute AIC value or number of parameters in the model (Marcoulides & Hershberger, 2014), with Δ AIC > 2 considered a significant difference (Burnham & Anderson, 2002). If model fit improved significantly, estimated coefficients of individual parameters were interpreted. If a model failed to converge it was simplified as appropriate. As therapeutic alliance ratings were measured at three time-points, only linear slopes were tested. Following the selection of the best fitting unconditional model, all hypothesized predictors were added simultaneously for the final model. If the inclusion of random slopes was previously shown to improve model fit (i.e., if growth trajectories differed significantly across therapeutic dyads), a model with cross-level interaction effects was also calculated. The full maximum likelihood estimator was used during model selection, however final models were calculated using the more stable restricted maximum likelihood estimator (Snijders & Bosker, 2011). Hypothesis testing was conducted with significance predefined as p < 0.05 and confidence intervals, as reliance on p-values is heavily debated in multilevel models (Bates, 2006). A Bonferroni correction to the p-value threshold for significance was applied in the final stage of analyses, to control for the increased possibility of Type I error associated with testing multiple outcome variables. After correcting for the five models calculated in the current study, the p-value indicating significance was set at p < .01. It should be noted that a number of familiar concepts from simple linear models do not translate to linear mixed models, including correction for multiple comparisons and R² (Luke, 2019). A number of pseudo-R² statistics for multilevel models have been proposed, however there has been little consensus as to the most accurate and appropriate calculation, definition and interpretation (LaHuis et al., 2014; Roberts et al., 2011). The most popular alternative attempts to summarize model goodness of fit with two statistics—the marginal R² which describes the variance explained by fixed effects, and the conditional R² which describes the proportion of variance explained by the overall model (Nakagawa & Schielzeth, 2013). These statistics are used in the current study, with the acknowledgement they are not analogous to traditional R² and should be interpreted and evaluated alongside other indicators of model fit (Jaeger et al., 2017; Nakagawa & Schielzeth, 2013).

Model assumptions were checked in R. Using Pearson correlation coefficients, the chi-square test and independent t-tests with Cohen’s d effect sizes, no multicollinearity was detected. Problematic outliers were defined as values that deviate markedly from other data-points in the sample, and were either errors (i.e., out of range values as opposed to “interesting” outliers that represent naturally occurring observations; Aguinis et al., 2013) or exerted undue influence over results. Cook’s distance plots were used to detect influential outliers using the “predictmeans” package in R (Luo et al., 2018), with a cut-off value of 0.50 (Van der Meer et al., 2010). Ten percent of all data were examined for accuracy in data entry, with no out of range values detected. While all therapeutic alliance subscales were somewhat skewed and contained outliers representing less common scores (e.g., a notably low level of client participation or high level of therapist negative attitude), no influential outliers were detected using a cut-off of 0.5 (there is debate over the appropriate cut-off value for Cook's D, with this cut-off representing a middle ground; Van der Meer et al., 2010). A degree of expected heteroscedasticity was noted in all models.

Results

Characteristics of the Sample

Relevant clinical characteristics for the total sample and by index illness are presented in Table I. Total Beck Depression Inventory scores were not recorded for two participants, and the exact duration of primary illness could not be coded for one participant in the MD study, so these data points are not included.

Table I. Clinical characteristics in the total sample and by diagnostic group.

			Diagnostic group
	Total sample (N = 212)		Depression (n = 100)		Bulimia nervosa (n = 58)		Binge eating disorder (n = 54)
	M/N	SD/%	M/n	SD/%	M/n	SD/%	M/n	SD/%
Primary illness onset	21.2	10.0	21.75	11.49	20.91	8.39	21.52	10.40
GAF	55.12	7.26	55.16	9.37	54.17	5.30	56.06	3.69
MADRS	18.91	9.84	24.31	6.48	15.33	9.58	12.76	10.05
Prior self-harm	58	27%	25	25%	20	35%	13	24%
Suicide attempt	46	22%	26	26%	13	22%	7	13%
BDI	26.28	9.66	26.28	9.66
Eating Disorder Inventory
Thinness drive	12.42	5.47			14.05	5.25	10.67	5.19
Bulimia	10.17	4.23			10.83	4.52	9.46	3.83
Dissatisfaction	19.40	7.80			19.02	8.27	19.81	7.32

Note. BDI: Beck Depression Inventory; GAF: Global Assessment of Functioning; MADRS: Montgomery-Asberg Depression Rating Scale.

Childhood trauma and perceived parental bonding among the sample

Table II displays childhood trauma and PBI subscale data for the total sample and by index illness. Thirteen participants did not complete the PBI for paternal caregivers, and one did not complete the maternal version of the same questionnaire, consistent with the absence of this caregiver in their lives. Two additional participants did not complete one subscale of the PBI for paternal caregivers.

Table II. Perceived parental bonding and childhood trauma in the total sample and by diagnostic group.

			Diagnostic group
	Total sample (N = 212)		Major depression (n = 100)		Bulimia nervosa (n = 58)		Binge eating disorder (n = 54)
	M/n	SD/%	M/n	SD/%	M/n	SD/%	M/n	SD/%
Maternal bonding
Care	20.96	9.85	21.28	9.14	21.98	10.80	19.28	10.02
Overprotection	15.09	8.62	15.31	8.61	15.00	9.43	14.78	7.89
Paternal bonding
Care	18.20	9.77	16.18^a	9.13	19.07^ab	10.95	20.75^b	8.95
Overprotection	13.03	7.47	13.21	7.85	12.43	7.38	13.35	6.99
Childhood trauma
Impersonal trauma	7	3%	5	5%	1	2%	1	2%
Psychological abuse	66	31%	29	29%	20	35%	17	32%
Physical abuse	50	23%	17	17%	12	20%	15	28%
Sexual abuse	63	30%	28	29%	20	34%	15	28%

Overall, the sample reported low paternal (M = 18.2, SD = 9.8) and maternal care scores (M = 21, SD = 9.9) and high paternal (M = 13, SD = 7.5) and maternal (M = 15.1, SD = 8.6) overprotection scores. The most common form of childhood trauma reported was psychological abuse (31% overall), followed by sexual abuse (30% overall: 4% non-contact, 12% contact and 14% attempted or completed intercourse), physical abuse (23% overall: 11% threatened, 6% occurred and 6% occurred and required medical attention) and witnessing or experiencing a disaster/accident/war (3%). Forty-seven percent of the sample did not report traumatic experiences before the age of 16, 31.5% reported one, 13.8% two, 6.6% three and 1.1% four or more, respectively. Severity, frequency and polyvictimisation data were not incorporated in modelling due to low base rates. Subscales of the PBI had weak to moderate correlations with one another, which did not approach the threshold considered to indicate problems with multicollinearity (> .70; Tabachnick & Fidell, 1996). Correlations among childhood trauma variables, and between parental bonding and childhood trauma variables, were similarly either moderate, weak or non-significant. Accidents and disasters were not associated with any other predictor variables. Participants who reported higher levels of care and lower levels of parental overprotection were less likely to report childhood abuse experiences.

Inter-rater reliability, inter-rater agreement and internal consistency

CVs, ICC estimates with 95% confidence intervals, coefficient alpha estimates with 95% confidence intervals and coefficient omega estimates were calculated for each subscale. Range restriction was present in some, especially those reflecting the contribution of therapist behaviour and demeanour to the overall therapeutic alliance. Scores were concentrated at the alliance-supportive end of these subscales (e.g., high scores on the therapist warmth subscale of the VPPS) and this lack of variability reduced inter-rater reliability and internal consistency estimates despite high inter-rater agreement (Hallgren, 2012; Trevethan, 2017). In the interests of detecting effects, models were calculated only for subscales that demonstrated at least moderate inter-rater reliability (the VTAS patient contribution and VPPS patient participation, patient distress, patient dependency and patient hostility subscales). The strongest correlations among therapeutic alliance subscales were between the VPPS patient participation and VTAS patient contribution subscales at all three time points (r₁= 0.68, r₂= 0.64, r₃ = 0.66), indicating these subscales captured somewhat similar aspects of patient behaviour in therapy. This association did not reach the threshold suggesting potential problems with multicollinearity, however (Tabachnick & Fidell, 1996).

Linear Mixed Models

Linear mixed models were calculated to characterize the mean and individual effects of time on client participation, dependency, hostility, and overall contribution to alliance, and to examine whether any observed variation was due to therapy type, primary diagnosis, childhood trauma or perceived parental bonding. Estimated coefficients and model fit indices for the final models are displayed in Table III.

Table III. Final linear mixed models for therapeutic alliance ratings in the total sample.

	VPPS patient participation	VPPS patient dependency	VPPS patient distress	VPPS patient hostility	VTAS patient contribution
	Estimates [CI]	Estimates [CI]	Estimates [CI]	Estimates [CI]	Estimates [CI]
Predictors
Intercept	4.26 [4.01, 4.50]**	2.07 [1.86,2.27]**	1.47 [1.20, 1.75]**	1.16 [1.06, 1.27]**	4.33 [4.10, 4.57]**
Phase of therapy	0.04 [0.006, 0.07]*	−0.005 [−0.02, 0.01]	−0.04 [−0.09, 0.006]	−0.008 [−0.08, 0.07]	0.06 [0.03, 0.10]**
Therapy (CBT)	0.07 [−0.06, 0.20]	0.003 [−0.11, 0.11]	−0.05 [−0.20, 0.10]	−0.05 [−0.10, 0.009]	0.07 [−0.06, 0.19]
Diagnosis (MD)	−0.08 [−0.21, 0.05]	−0.39 [−0.50, −0.28]**	0.36 [0.21, 0.51]**	−0.03 [−0.09, 0.22]	−0.18 [−0.31, −0.05]**
Childhood trauma
Impersonal trauma	−0.28 [−0.65, 0.09]	−0.03 [−0.34, 0.28]	−0.12 [−0.53, 0.30]	0.03 [−0.14, 0.19]	−0.05 [−0.41, 0.31]
Physical abuse	−0.06 [−0.23, 0.10]	0.03 [−0.11, 0.17]	−0.02 [−0.20, 0.16]	−0.002 [−0.07, 0.07]	−0.10 [−0.26, 0.05]
Sexual abuse	0.07 [−0.07, 0.22]	−0.04 [−0.16, 0.08]	−0.04 [−0.20, 0.12]	0.03 [−0.03, 0.09]	0.04 [−0.10, 0.18]
Psychological abuse	0.16 [−0.009, 0.33]	0.04 [−0.10, 0.18]	0.12 [−0.06, 0.31]	0.01 [−0.06, 0.09]	0.12 [−0.04, 0.28]
PBI
Maternal care	−0.002 [−0.01, 0.006]	−0.004 [−0.01, 0.002]	0.003 [−0.005, 0.01]	−0.0003 [−0.004, 0.003]	0.002 [−0.006, 0.01]
Maternal OP	−0.002 [−0.01, 0.008]	−0.004 [−0.01, 0.004]	0.001 [−0.01, 0.01]	−0.0008 [−0.005, 0.004]	0.003 [−0.006, 0.01]
Paternal care	0.009 [0.001, 0.02]*	−0.004 [−0.01, 0.002]	−0.003 [−0.01, 0.005]	0.0008 [−0.002, 0.004]	0.004 [−0.003, 0.01]
Paternal OP	0.007 [−0.003, 0.02]	−0.0003 [−0.009, 0.008]	0.001 [−0.01, 0.01]	−0.0001 [−0.004, 0.004]	0.003 [−0.006, 0.01]
Interactions
Therapy (CBT) x phase				−0.004 [−0.04, 0.04]
Diagnosis (MD) x phase				0.03 [−0.01, 0.07]
Impersonal trauma x phase				−0.02 [−0.13, 0.09]
Physical abuse x phase				0.01 [−0.04, 0.06]
Sexual abuse x phase				−0.05 [−0.09, −0.004]*
Psychological abuse x phase				−0.003 [−0.05, 0.004]
Maternal care x phase				−0.002 [−0.004, 0.0008]
Maternal OP x phase				−0.001 [−0.004, 0.002]
Paternal care x phase				−0.0007 [−0.003, 0.002]
Paternal OP x phase				0.001 [−0.001, 0.004]
Random effects
Intercept standard deviation	0.37 [0.32, 0.42]	0.21 [0.07, 0.56]	0.38 [0.33, 0.45]	0.11 [0.08, 0.15]	0.35 [0.30, 0.40]
Slope standard deviation				0.06 [0.28, 0.11]
Intercept-slope correlation				0.16 [−0.9, 0.99]
Residual	0.03	0.44	0.38	0.03	0.32
Model Comparison
R^2 marginal	0.07	0.15	0.12	0.04	0.07
R^2 conditional	0.61	0.30	0.57	0.52	0.57
Observations (N)	504 (168)	504 (168)	504 (168)	504(168)	504 (168)

Note. CI: Confidence intervals; CBT: cognitive behavioural therapy; MD: Major depression; OP: Overprotection, *p < 0.05, **p < 0.01.

Subscale results

VPPS Patient Participation Subscale. The mean participation score was 4.26 (se = 0.03). The ICC from a random-intercept model was 0.59, indicating that 59% of the total variation in initial scores was between therapeutic dyads, with the remainder of variability occurring within therapeutic dyads across time. A random linear slope was not significantly better than a fixed linear slope, indicating data were best modelled by an average growth trajectory that increased by 0.04 at each time point. The default independent error structure was retained. In the final model, the expected value at time one was 4.26 (se = 0.12) and none of the hypothesized predictor variables accounted for variation in the intercept or slope.

VPPS Patient Dependency Subscale. The mean dependency score was 1.48 (se = 0.02) and the ICC from a random-intercept model was 0.61. The inclusion of fixed or random linear slopes did not improve the model, suggesting there were no significant changes in VPPS patient dependency scores over time. A fixed linear slope was retained for theoretical purposes, however. An unstructured error-covariance structure produced the smallest AIC, indicating residuals differed between therapeutic dyads across observations. In the final model, the expected value at time one was 2.07 (se = 0.10). The estimated primary diagnosis coefficient (β = −0.39, se = 0.06) indicated that participants diagnosed with depression had an initial dependency score 0.39 below those diagnosed with bulimia nervosa or binge eating disorder.

VPPS Patient Distress Subscale. The mean distress score was 1.94 (se = 0.03) and the ICC from a random-intercept model was 0.45. A model with random linear slopes was not significantly better than one with a fixed linear slope, indicating data were best modelled by an average growth trajectory that decreased by 0.06 at each time point. A model with an unstructured error-covariance structure produced the smallest AIC. In the final model, the expected value at time one was 1.48 (se = 0.13). Primary diagnosis was significantly related to differences in initial scores, t(157) = 4.86, p < 0.01. The estimated coefficient (β = .36, se = 0.07) indicated that participants with a primary diagnosis of mood disorder had an initial distress score 0.36 above those diagnosed with bulimia nervosa or binge eating disorder.

VPPS Patient Hostility Subscale. The mean hostility score was 1.11 (se = 0.009) and the ICC from a random-intercept model was 0.42. A model with a fixed linear slope indicated an increase in hostility of 0.02 at each time point. A model with random linear slopes produced a significantly lower AIC than a model with a fixed slope, indicating data were better described by a model allowing the initial level and development of hostility scores over time to vary across therapeutic dyads. The standard deviation of random slopes was 0.05. There was no significant correlation between random intercepts and random slopes, indicating there was no association between a client’s initial level of hostility in therapy and time. In the final two models, where main effects and interaction effects were added sequentially with a random slope retained, none of the hypothesized predictor variables accounted for variation in the intercept or slope.

VTAS Patient Contribution Subscale. The mean patient contribution score was 4.19 (se = 0.02) and the ICC from a random-intercept model was 0.54. A random linear slope was not significantly better than a fixed linear slope, indicating data were best modelled by an average growth trajectory that increased by 0.06 at each time point. The best-fitting error covariance structure was independent. In the final model, the expected value at time one was 4.33 (se = 0.12). Primary diagnosis was significantly related to differences in initial patient contribution scores, t(157) = −2.80, p < 0.01. The estimated coefficient (β = −0.18, se = 0.06) indicated that participants with a primary diagnosis of major depression had an initial score 0.18 below those diagnosed with bulimia nervosa or binge eating disorder.

Discussion

The current study examined therapeutic alliance in the combined MD and BE sample, including assessing whether childhood and clinical factors influenced alliance development. Mean alliance scores in early therapy, as measured by both the VPPS and VTAS, indicated that participants formed strong alliances with their therapists overall. More specifically, high levels of participation in therapeutic tasks and other alliance-promoting behaviours, moderate distress and only slight hostility towards and dependence on therapists were recorded by raters. Models indicated significant variability among participants’ initial scores on all subscales, but similar growth trajectories over time in all but the patient hostility subscale. Client-driven aspects of the alliance improved on average as therapy progressed, with the exception of patient dependency scores which did not change significantly and patient hostility scores which increased slightly.

The hypothesized associations between poor parental bonds and childhood trauma with weaker therapeutic alliances were not supported. It is possible this may be due in part to necessary though unfortunately harsh adjustments for multiple comparisons. While the use of the Bonferroni correction reduced the risk of Type I error in the study, it also increased the risk of Type II error (Bates, 2006; Feise, 2002). For example, the Bonferroni correction assumes independence among tests, an assumption violated by the weak to moderate correlations found among therapeutic alliance subscale scores at all three time points in therapy. Applying a less conservative statistical adjustment to the p-value threshold for significance was not possible given the complexity of multilevel models, and applying adjustments to individual p-values in all five models posed further danger of over-correction, since multiple comparisons within the linear mixed models were already accounted for by model calculation and selection procedures (Gelman et al., 2012). The hypothesized predictor variables thus continue to warrant further exploration.

It was additionally hypothesized that any deleterious influences on therapeutic alliance strength in early therapy accounted for by predictor variables would lessen over time, with this effect being more pronounced in schema therapy compared with CBT. Differences in therapeutic alliance strength across therapy modalities have typically been explored via meta-analysis (Flückiger et al., 2018), however the design of the current study allowed for the direct comparison of two psychotherapies in the MD sample and three in the BE sample. Although CBT conceptualizes the alliance as a supportive framework rather than a direct intervention (Beck et al., 1979; Raykos et al., 2014), no significant differences were found between CBT and schema therapy in alliance formation or alliance growth trajectory for any of the measured subscales. As schema therapy is a form of CBT, shared elements across the two therapies could account for this lack of difference in alliance scores. Future studies directly comparing dissimilar therapies, such as CBT and psychodynamic or interpersonal psychotherapy, would provide greater clarity about potential differences in therapeutic alliance development related to different schools of therapy.

Primary diagnosis was a significant predictor of difference in several initial alliance subscale scores, with a diagnosis of depression being associated with higher distress, lower dependency and a lesser client-driven contribution to a strong therapeutic alliance. Binge eating disorder and bulimia nervosa were in turn associated with lower distress, higher dependency and greater client-driven contribution to a strong therapeutic alliance. Findings suggest the possibility of unique challenges in building the alliance with different clinical populations, potentially due to associated interpersonal patterns or personality traits, such as withdrawal from social interactions when depressed (Abbate-Daga et al., 2013; Black et al., 2013). Early awareness of potential difficulties in alliance formation may assist therapists to better meet their clients’ alliance needs and thus maximize treatment outcomes, by selecting therapies that target alliance difficulties as part of treatment (Zilcha-Mano et al., 2021) or by adjusting therapist relational style or task selection with diagnosis-specific alliance challenges in mind. For instance, by responding early to signs of unhelpful dependency when working with clients diagnosed with binge eating disorder or bulimia nervosa.

All final models accounted for a small proportion of variance in therapeutic alliance scores, as indicated by marginal and conditional R² values. The small effect sizes observed are likely linked to several limitations of the current study. Challenges in accurately measuring childhood trauma and parental bonds, including inconsistent or under reporting and subjective interpretation of events, are well documented and apply here as in any research of the topic (e.g., Frissa et al., 2016). The practical demands of therapeutic alliance data collection via independent observer resulted in the rating of only three time-points, limiting the trends tested in the models to fixed and random linear slopes when it is possible alliance development over therapy in the study would be better represented by cubic, plateau or other higher order polynomial trends.. Future research taking a similar approach to studying therapeutic alliance should endeavour to use larger sample sizes and assess therapeutic alliance across a greater number of time-points, perhaps necessitating the use of less time-intensive therapeutic alliance measures (such as self- or therapist-report) despite their potential biases (Cecero et al., 2001; Elvins & Green, 2008).

Lack of variability in therapeutic alliance scores also limited the number of subscales tested, reduced the likelihood of detecting significant effects and prevented examination of individuals’ state-like alliance changes for most subscales. Subscales reflecting therapist behaviour and demeanour in particular were highly restricted, leading to poor inter-rater reliability estimates (as these assess the ability to distinguish cases as well as agreement between independent raters; Hallgren, 2012; Quarfoot & Levine, 2016). Clients displayed greater variation in alliance scores, although relatively high levels of alliance-promoting and low levels of alliance-undermining behaviours were still observed overall. Range restriction was likely linked to sample characteristics, such as the level of client motivation required to participate in a year-long RCT. Further, that therapists were well trained, experienced and closely supervised likely contributed to their consistent demonstrations of alliance-supportive behaviour, such as extremely low levels of negative attitudes and high levels of warmth towards clients. Future research would likely benefit from building variation representative of “real world” therapy into the study design. This could be achieved by selecting therapists early in their career, as it can be assumed that some novice therapists will have more difficulty managing alliance ruptures compared with their more experienced counterparts (Bilican & Soygut, 2015; Callahan & Hynan, 2005; Cartwright et al., 2014; Kline et al., 2019; Swift & Greenberg, 2012). Additionally, clients with a degree of ambivalence about engaging in treatment, such as those in mandated treatment settings, may produce more fruitful data for answering the question of how client characteristics and early experiences impact therapeutic alliance development. Another possible explanation for the blunted therapeutic alliance scores observed in the current study is the use of observer-rated therapeutic alliance measures, as these are inherently unable to capture internal aspects of the alliance (Elvins & Green, 2008). The use of audio-recordings as the basis for ratings, as opposed to video-recordings, may have also limited the ability of raters to detect nuances in alliance strength. As well as reducing the time commitment involved in making independent ratings, self-report measures may therefore produce more useful data.

A major strength of the current research was the assessment of therapeutic alliance over three time-points and the use of linear mixed modelling to examine the contribution of trait and state-like differences in therapeutic alliance (Zilcha-Mano & Fisher, 2022). The effect of predictor variables and state-like changes in alliance scores were unable to be adequately explored due to low variation in alliance data, however differences in trait-like alliance related to diagnosis nonetheless provide support for the premise that certain client subgroups may begin therapy with a lesser chance of forming an immediately strong therapeutic alliance, and will benefit from therapists who are aware of this and select treatment with an individual’s alliance needs in mind. Further research that captures the dynamic and individualized nature of therapeutic alliance development is needed to better understand the associations between childhood experiences, clinical characteristics and therapeutic alliance over different therapy phases.

Ethics Statement

The Depression and Binge Eating studies both obtained ethical approval from the relevant bodies and were lodged with the Australian and New Zealand Clinical Trials Registry (depression study registry number: ACTRN12605000723684; binge eating study registry number: ACTRN12605000721606). Additionally, this study has sought and received approval from the University of Canterbury. Human Research Ethics Committee. All participants provided informed written consent for their participation.

Acknowledgements

We acknowledge study coordinators Helen Kleindienst, Caroline Bray, Sarah Rowe, Andrea Bartram and Julia Martin.

Disclosure Statement

No potential conflict of interest was reported by the author(s).

Data Availability Statement

The data included in this manuscript have not been made publicly available for ethical and privacy reasons.

Supplemental Data

Supplemental data for this article can be accessed at https://doi.org/10.1080/10503307.2023.2191800.

Additional information

Funding

This research was supported by funding from a University of Canterbury Doctoral Scholarship, a Lottery Health Grant, and a Programme Grant from the Health Research Council of New Zealand [grant number HRC04/292B].