Overview

This document applies the same preregistered exclusion criteria used in the Study 1 full-data analysis to the new-design pilot data. The DVs themselves differ from Study 1: rather than separate retrospective (“past week”) and prospective/forecasted item pairs per construct, this design collects a single forecast item per construct, anchored relative to a normal week of companion use (e.g., “Much less lonely… No different… Much more lonely”). There are no reverse-coded items in this design (all anchors are already oriented so that responses can be used directly), so no reverse-scoring step is needed before constructs are formed.

# install.packages(c("qualtRics", "dplyr", "tidyr", "purrr", "broom", "ggplot2"))
library(qualtRics)
library(dplyr)
library(tidyr)
library(purrr)
library(broom)
library(ggplot2)

0. Import data

Update this path to point to your pilot export (standard Qualtrics CSV export — legacy format with the three header rows: field name, question text, ImportId JSON).

pilot <- read_survey("~/Google drive/My Drive/YEAR 2/PROJECTS/SANDRA/AI Companionship/Cold Turkey Study/Expectations/Data/newdesign_pilot.csv")

cat("Pilot N (raw):", nrow(pilot), "\n")
## Pilot N (raw): 88

1. Remove preview / test responses

Same as the Study 1 script — preview responses (from testing the survey link) are dropped before applying substantive exclusion criteria.

pilot <- pilot %>% filter(DistributionChannel != "preview")

cat("N after removing preview responses:", nrow(pilot), "\n")
## N after removing preview responses: 85

2. Preregistered exclusion criteria (same protocol as Study 1, Section 6)

Participants are excluded if they do not pass a standard attention-check, Qualtrics’ bot detection software, a hidden-text bot detection item, Qualtrics’ duplicate-response detection, and/or if they complete the survey in less than 1 minute. Operationalized identically to the Study 1 script:

Respondents screened out as ineligible before reaching the substantive survey blocks (i.e., never reached consent, and therefore missing all well-being data) are dropped as part of this step as well.

clean_pilot <- pilot %>%
  mutate(
    reached_survey  = !is.na(consent),
    is_duplicate    = ifelse(is.na(Q_DuplicateRespondent), FALSE, Q_DuplicateRespondent),
    fails_recaptcha = Q_RecaptchaScore < 0.5
  ) %>%
  filter(
    reached_survey,
    `Duration (in seconds)` >= 60,
    attn == 2,
    !is.na(attn_bots),
    !is_duplicate,
    !fails_recaptcha
  )

cat("Final analytic N after exclusions:", nrow(clean_pilot), "\n")
## Final analytic N after exclusions: 27

3. Sample characteristics

clean_pilot %>% count(gender, sort = TRUE) %>% knitr::kable()
gender n
2 13
1 12
3 1
4 1 
clean_pilot %>% count(most_platform, sort = TRUE) %>% knitr::kable()
most_platform n
5 19
1 2
10 2
11 2
12 1
14 1 
clean_pilot %>% count(rel_status, sort = TRUE) %>% knitr::kable()
rel_status n
1 10
3 8
2 5
4 2
5 2

Note: these come out as numeric codes here since only the numeric export was provided for this pilot. If you export a text-labeled version of this same file (choice text instead of numeric values), it can be left-joined in on ResponseId the same way the Study 1 script merges txt in, to get readable labels for these categorical/demographic variables.

4. Construct scores

Each construct here is a single forecast item, not an average of two items, and none require reverse-scoring. Higher values always mean “more [construct label]” as literally anchored in the item (e.g., higher lonely = more lonely; higher selfesteem = better about yourself; higher connect = more connected) — this is not a consistent better/worse direction across constructs, so keep that in mind for any composite index built later (the same sign-flip logic Study 1 used in Section 10d, via a higher_is_better lookup, would apply here too).

construct_items <- c(
  loneliness        = "lonely",
  depression         = "depressed",
  irq                = "irq",
  happiness          = "happiness",
  life_satisfaction  = "swl",
  stress             = "stress",
  self_esteem        = "selfesteem",
  anxiety            = "anxiety",
  connectedness      = "connect"
)

for (construct in names(construct_items)) {
  item_col <- construct_items[[construct]]
  clean_pilot[[paste0(construct, "_forecast")]] <- clean_pilot[[item_col]]
}

higher_is_better <- c(
  loneliness = FALSE, depression = FALSE, stress = FALSE, anxiety = FALSE,
  happiness = TRUE, life_satisfaction = TRUE, self_esteem = TRUE,
  connectedness = TRUE, irq = TRUE
)

clean_pilot %>%
  select(ResponseId, ends_with("_forecast")) %>%
  summarise(across(ends_with("_forecast"),
                    list(mean = ~mean(.x, na.rm = TRUE),
                         sd   = ~sd(.x, na.rm = TRUE),
                         n    = ~sum(!is.na(.x))))) %>%
  pivot_longer(everything(),
               names_to = c("construct", "stat"),
               names_pattern = "(.*)_forecast_(mean|sd|n)") %>%
  pivot_wider(names_from = stat, values_from = value) %>%
  knitr::kable(digits = 2)
construct mean sd n
loneliness 5.11 1.19 27
depression 4.63 1.15 27
irq 4.96 1.16 27
happiness 3.04 0.98 27
life_satisfaction 3.52 1.19 27
stress 5.30 0.99 27
self_esteem 3.52 0.98 27
anxiety 5.11 1.01 27
connectedness 4.56 1.01 27

4a. Descriptives: difficulty, missing, and interaction-change items

Before folding these into the secondary analysis, a quick look at their distributions on their own: f_difficult (1 = Very easy to 5 = Very difficult), f_miss (1 = Not at all to 5 = A great deal), and the three interaction-change items — f_interaction (in-person), f_online, and f_alone — each 1 = Much less to 5 = Much more than normal, with 3 = About the same.

exploratory_dvs <- c("f_difficult", "f_miss", "f_interaction", "f_online", "f_alone")

exploratory_dv_summary <- map_dfr(exploratory_dvs, function(var) {
  x <- clean_pilot[[var]]
  tibble(
    variable = var,
    n        = sum(!is.na(x)),
    mean     = mean(x, na.rm = TRUE),
    sd       = sd(x, na.rm = TRUE),
    median   = median(x, na.rm = TRUE),
    min      = min(x, na.rm = TRUE),
    max      = max(x, na.rm = TRUE)
  )
})

knitr::kable(exploratory_dv_summary, digits = 2)
variable n mean sd median min max
f_difficult 27 3.33 1.24 4 1 5
f_miss 27 3.44 1.09 3 1 5
f_interaction 27 3.44 0.80 3 2 5
f_online 27 3.44 0.80 3 2 5
f_alone 27 3.19 0.74 3 1 4

Response distributions (count and % at each scale point):

exploratory_dv_dist <- map_dfr(exploratory_dvs, function(var) {
  clean_pilot %>%
    filter(!is.na(.data[[var]])) %>%
    count(value = .data[[var]]) %>%
    mutate(variable = var, pct = round(100 * n / sum(n), 1)) %>%
    select(variable, value, n, pct)
})

knitr::kable(exploratory_dv_dist, digits = 1)
variable value n pct
f_difficult 1 3 11.1
f_difficult 2 4 14.8
f_difficult 3 5 18.5
f_difficult 4 11 40.7
f_difficult 5 4 14.8
f_miss 1 1 3.7
f_miss 2 4 14.8
f_miss 3 9 33.3
f_miss 4 8 29.6
f_miss 5 5 18.5
f_interaction 2 3 11.1
f_interaction 3 11 40.7
f_interaction 4 11 40.7
f_interaction 5 2 7.4
f_online 2 3 11.1
f_online 3 11 40.7
f_online 4 11 40.7
f_online 5 2 7.4
f_alone 1 1 3.7
f_alone 2 2 7.4
f_alone 3 15 55.6
f_alone 4 9 33.3 
ggplot(exploratory_dv_dist, aes(x = factor(value), y = pct)) +
  geom_col(fill = "steelblue") +
  facet_wrap(~variable, ncol = 2) +
  labs(x = "Scale point (1–5)", y = "% of respondents",
       title = "Distributions: Difficulty, Missing, and Interaction-Change Items") +
  theme_minimal()

5. Primary analysis: does the average forecast differ from “no different”?

This design has no retrospective baseline to difference against (unlike Study 1’s paired retro-vs-forecast test), so the analogous primary test here is a one-sample t-test comparing each construct’s forecast score against the scale’s neutral midpoint — 4 on the 1–7 scale, i.e., “No different” — which functions as the null of “no expected change.”

midpoint <- 4

primary_results <- map_dfr(names(construct_items), function(construct) {
  forecast <- clean_pilot[[paste0(construct, "_forecast")]]
  dev      <- forecast - midpoint

  test <- t.test(forecast, mu = midpoint)

  tibble(
    construct              = construct,
    n                      = sum(!is.na(forecast)),
    mean_forecast          = mean(forecast, na.rm = TRUE),
    mean_dev_from_midpoint = mean(dev, na.rm = TRUE),
    t                      = unname(test$statistic),
    df                     = unname(test$parameter),
    p                      = test$p.value,
    dz                     = mean(dev, na.rm = TRUE) / sd(dev, na.rm = TRUE)
  )
}) %>%
  mutate(q_BH = p.adjust(p, method = "BH"))

knitr::kable(primary_results, digits = 3)
construct n mean_forecast mean_dev_from_midpoint t df p dz q_BH
loneliness 27 5.111 1.111 4.862 26 0.000 0.936 0.000
depression 27 4.630 0.630 2.849 26 0.008 0.548 0.011
irq 27 4.963 0.963 4.315 26 0.000 0.830 0.000
happiness 27 3.037 -0.963 -5.107 26 0.000 -0.983 0.000
life_satisfaction 27 3.519 -0.481 -2.105 26 0.045 -0.405 0.045
stress 27 5.296 1.296 6.784 26 0.000 1.306 0.000
self_esteem 27 3.519 -0.481 -2.565 26 0.016 -0.494 0.019
anxiety 27 5.111 1.111 5.701 26 0.000 1.097 0.000
connectedness 27 4.556 0.556 2.850 26 0.008 0.549 0.011

mean_dev_from_midpoint is signed exactly as the item is coded — e.g., a positive value for loneliness means the average forecast leans toward “more lonely,” while a positive value for happiness means the average forecast leans toward “more happy.” Remember higher_is_better (defined in Section 4) when interpreting whether a given sign is good or bad news for that construct — it is not consistent across rows of this table.

5a. % of participants predicting improvement, no change, or worsening

Using the same benefit/harm direction convention as Study 1’s Section 10d, via higher_is_better:

expectation_categories <- map_dfr(names(construct_items), function(construct) {
  dev         <- clean_pilot[[paste0(construct, "_forecast")]] - midpoint
  better_dir  <- higher_is_better[[construct]]

  category <- case_when(
    is.na(dev) ~ NA_character_,
    dev == 0 ~ "No change",
    (dev > 0 & better_dir) | (dev < 0 & !better_dir) ~ "Predicted improvement",
    TRUE ~ "Predicted worsening"
  )

  tibble(construct = construct, category = category)
})

expectation_summary <- expectation_categories %>%
  filter(!is.na(category)) %>%
  count(construct, category) %>%
  group_by(construct) %>%
  mutate(pct = round(100 * n / sum(n), 1)) %>%
  ungroup()

knitr::kable(expectation_summary, digits = 1)
construct category n pct
anxiety No change 6 22.2
anxiety Predicted improvement 1 3.7
anxiety Predicted worsening 20 74.1
connectedness No change 9 33.3
connectedness Predicted improvement 14 51.9
connectedness Predicted worsening 4 14.8
depression No change 10 37.0
depression Predicted improvement 2 7.4
depression Predicted worsening 15 55.6
happiness No change 11 40.7
happiness Predicted worsening 16 59.3
irq No change 6 22.2
irq Predicted improvement 19 70.4
irq Predicted worsening 2 7.4
life_satisfaction No change 12 44.4
life_satisfaction Predicted improvement 3 11.1
life_satisfaction Predicted worsening 12 44.4
loneliness No change 4 14.8
loneliness Predicted improvement 1 3.7
loneliness Predicted worsening 22 81.5
self_esteem No change 13 48.1
self_esteem Predicted improvement 3 11.1
self_esteem Predicted worsening 11 40.7
stress No change 6 22.2
stress Predicted worsening 21 77.8 
library(ggplot2)
ggplot(expectation_summary, aes(x = construct, y = pct, fill = category)) +
  geom_col(position = "stack") +
  coord_flip() +
  labs(x = NULL, y = "% of participants", fill = NULL,
       title = "Distribution of Forecasted Change by Construct") +
  theme_minimal()

6. Secondary analysis: difficulty, missing, and interaction change as predictors

Three additional single-item DVs were collected alongside the construct forecasts (see the “Exploratory DVs” block): how difficult the participant thinks the week would be (f_difficult, 1 = Very easy to 5 = Very difficult), how much they think they’d miss their companion (f_miss, 1 = Not at all to 5 = A great deal), and three items on anticipated changes to time spent interacting with others in-person (f_interaction), online (f_online), and alone (f_alone) — each on a 1 = Much less to 5 = Much more scale, relative to normal.

One regression per construct per predictor, using the raw forecast deviation from the midpoint as the outcome (same sign convention as Section 5 — not re-oriented for harm/benefit here, so interpret sign per-construct as above). All unadjusted, exploratory p-values.

secondary_predictors <- c("f_difficult", "f_miss", "f_interaction", "f_online", "f_alone")

secondary_results <- map_dfr(names(construct_items), function(construct) {
  dev <- clean_pilot[[paste0(construct, "_forecast")]] - midpoint

  map_dfr(secondary_predictors, function(predictor) {
    pred <- clean_pilot[[predictor]]
    model <- lm(dev ~ pred)
    coefs <- broom::tidy(model)
    slope_row <- coefs[coefs$term == "pred", ]

    tibble(
      construct = construct,
      predictor = predictor,
      n         = sum(!is.na(dev) & !is.na(pred)),
      estimate  = slope_row$estimate,
      se        = slope_row$std.error,
      p         = slope_row$p.value
    )
  })
}) %>%
  group_by(predictor) %>%
  mutate(q_BH = p.adjust(p, method = "BH")) %>%
  ungroup()

knitr::kable(secondary_results, digits = 3)
construct predictor n estimate se p q_BH
loneliness f_difficult 27 0.300 0.182 0.111 0.201
loneliness f_miss 27 0.674 0.172 0.001 0.006
loneliness f_interaction 27 0.520 0.278 0.073 0.219
loneliness f_online 27 0.520 0.278 0.073 0.256
loneliness f_alone 27 0.387 0.313 0.229 0.682
depression f_difficult 27 0.208 0.180 0.259 0.333
depression f_miss 27 0.471 0.189 0.020 0.037
depression f_interaction 27 0.327 0.279 0.253 0.380
depression f_online 27 0.447 0.273 0.114 0.256
depression f_alone 27 0.203 0.310 0.519 0.682
irq f_difficult 27 0.083 0.186 0.658 0.658
irq f_miss 27 -0.083 0.213 0.699 0.699
irq f_interaction 27 0.567 0.267 0.044 0.196
irq f_online 27 0.267 0.285 0.358 0.537
irq f_alone 27 -0.271 0.311 0.391 0.682
happiness f_difficult 27 -0.433 0.132 0.003 0.014
happiness f_miss 27 -0.471 0.154 0.005 0.016
happiness f_interaction 27 -0.027 0.245 0.914 0.914
happiness f_online 27 -0.147 0.243 0.552 0.709
happiness f_alone 27 0.342 0.257 0.196 0.682
life_satisfaction f_difficult 27 -0.217 0.187 0.257 0.333
life_satisfaction f_miss 27 -0.203 0.215 0.355 0.456
life_satisfaction f_interaction 27 0.407 0.286 0.167 0.300
life_satisfaction f_online 27 -0.073 0.297 0.807 0.908
life_satisfaction f_alone 27 0.242 0.319 0.456 0.682
stress f_difficult 27 0.483 0.128 0.001 0.008
stress f_miss 27 0.536 0.148 0.001 0.006
stress f_interaction 27 0.147 0.246 0.557 0.716
stress f_online 27 0.327 0.239 0.184 0.332
stress f_alone 27 -0.034 0.270 0.900 0.900
self_esteem f_difficult 27 -0.267 0.148 0.084 0.188
self_esteem f_miss 27 -0.072 0.179 0.689 0.699
self_esteem f_interaction 27 0.347 0.234 0.150 0.300
self_esteem f_online 27 -0.013 0.244 0.957 0.957
self_esteem f_alone 27 0.100 0.264 0.708 0.797
anxiety f_difficult 27 0.300 0.152 0.059 0.178
anxiety f_miss 27 0.413 0.167 0.021 0.037
anxiety f_interaction 27 0.100 0.252 0.695 0.782
anxiety f_online 27 0.400 0.240 0.108 0.256
anxiety f_alone 27 0.174 0.273 0.531 0.682
connectedness f_difficult 27 0.125 0.161 0.446 0.502
connectedness f_miss 27 0.337 0.174 0.064 0.096
connectedness f_interaction 27 0.620 0.221 0.009 0.085
connectedness f_online 27 0.560 0.227 0.021 0.187
connectedness f_alone 27 0.300 0.269 0.275 0.682

6a. Summary of significant predictors (p < .05, unadjusted)

sig_secondary <- secondary_results %>%
  filter(p < .05) %>%
  arrange(predictor, p)

knitr::kable(sig_secondary, digits = 3)
construct predictor n estimate se p q_BH
stress f_difficult 27 0.483 0.128 0.001 0.008
happiness f_difficult 27 -0.433 0.132 0.003 0.014
connectedness f_interaction 27 0.620 0.221 0.009 0.085
irq f_interaction 27 0.567 0.267 0.044 0.196
loneliness f_miss 27 0.674 0.172 0.001 0.006
stress f_miss 27 0.536 0.148 0.001 0.006
happiness f_miss 27 -0.471 0.154 0.005 0.016
depression f_miss 27 0.471 0.189 0.020 0.037
anxiety f_miss 27 0.413 0.167 0.021 0.037
connectedness f_online 27 0.560 0.227 0.021 0.187 
if (nrow(sig_secondary) > 0) {
  cat("\n")
  for (i in seq_len(nrow(sig_secondary))) {
    row <- sig_secondary[i, ]
    direction <- ifelse(row$estimate > 0, "higher", "lower")
    cat(sprintf(
      "- %s predicts %s forecasted deviation in %s (b = %.3f, p = %.3f, n = %d)\n",
      row$predictor, direction, row$construct, row$estimate, row$p, row$n
    ))
  }
} else {
  cat("No secondary predictor effects reached p < .05.\n")
}
## 
## - f_difficult predicts higher forecasted deviation in stress (b = 0.483, p = 0.001, n = 27)
## - f_difficult predicts lower forecasted deviation in happiness (b = -0.433, p = 0.003, n = 27)
## - f_interaction predicts higher forecasted deviation in connectedness (b = 0.620, p = 0.009, n = 27)
## - f_interaction predicts higher forecasted deviation in irq (b = 0.567, p = 0.044, n = 27)
## - f_miss predicts higher forecasted deviation in loneliness (b = 0.674, p = 0.001, n = 27)
## - f_miss predicts higher forecasted deviation in stress (b = 0.536, p = 0.001, n = 27)
## - f_miss predicts lower forecasted deviation in happiness (b = -0.471, p = 0.005, n = 27)
## - f_miss predicts higher forecasted deviation in depression (b = 0.471, p = 0.020, n = 27)
## - f_miss predicts higher forecasted deviation in anxiety (b = 0.413, p = 0.021, n = 27)
## - f_online predicts higher forecasted deviation in connectedness (b = 0.560, p = 0.021, n = 27)