3. Merge numeric + text data
The following categorical/demographic variables are replaced with
their text-labeled values (from txt), since numeric codes
aren’t meaningful on their own for these fields (e.g., platform
selections, relationship framing, demographics). All scale items used
for construct scoring, and all variables used in the exclusion criteria,
are left as numeric (from num).
categorical_vars <- c(
"ai_purposes", "companion_frequency", "platforms", "most_platform",
"relation", "gender", "race", "ses", "employment", "edu", "overall_poli",
"rel_status", "have_child", "live_alone", "prior_dis_occurrence",
"companion_def_use", "unsure_followup", "def_frequency"
)
txt_subset <- txt %>% select(ResponseId, all_of(categorical_vars))
clean_full <- num %>%
select(-all_of(categorical_vars)) %>%
left_join(txt_subset, by = "ResponseId")
4. Preregistered exclusion criteria (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 as:
attn — instructional manipulation
check (“Please select ‘Somewhat disagree’”); correct response = 2
Q_RecaptchaScore — Qualtrics’ built-in
bot detection (reCAPTCHA v3); excluded if score < 0.5
Q_DuplicateRespondent — Qualtrics’
built-in duplicate-response detection; excluded if flagged
TRUE
Duration (in seconds) < 60 —
excluded as too fast to answer thoughtfully
Respondents screened out as ineligible before reaching the
substantive survey blocks (and therefore missing all well-being data)
are dropped as part of this step as well.
clean_full <- clean_full %>%
mutate(
attn_bots_clean = tolower(trimws(attn_bots)),
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_full), "\n")
## Final analytic N after exclusions: 214
6. Construct scores
Reverse-score the two “(R)” items, then average the two items per
construct to form a construct-level score (range 1-7), separately for
the Retrospective (Block A, r_) and Forecasted (Block B,
f_) item sets.
reverse_7pt <- function(x) 8 - x
clean_full <- clean_full %>%
mutate(
r_rse_1R = reverse_7pt(r_rse_1R),
f_rse_1R = reverse_7pt(f_rse_1R),
r_connect_1R = reverse_7pt(r_connect_1R),
f_connect_1R = reverse_7pt(f_connect_1R)
)
construct_items <- list(
loneliness = list(r = c("r_loneliness_1", "r_loneliness_2"),
f = c("f_loneliness_1", "f_loneliness_2")),
depression = list(r = c("r_depression_1", "r_depression_2"),
f = c("f_depression_1", "f_depression_2")),
anxiety = list(r = c("r_anxiety_1", "r_anxiety_2"),
f = c("f_anxiety_1", "f_anxiety_2")),
stress = list(r = c("r_stress_1", "r_stress_2"),
f = c("f_stress_1", "f_stress_2")),
happiness = list(r = c("r_happiness_1", "r_happiness_2"),
f = c("f_happiness_1", "f_happiness_2")),
life_satisfaction = list(r = c("r_swl_1", "r_swl_2"),
f = c("f_swl_1", "f_swl_2")),
self_esteem = list(r = c("r_rse_1R", "r_rse_2"),
f = c("f_rse_1R", "f_rse_2")),
irq = list(r = c("r_IRQ_1", "r_IRQ_2"),
f = c("f_IRQ_1", "f_IRQ_2")),
connectedness = list(r = c("r_connect_1R", "r_connect_2"),
f = c("f_connect_1R", "f_connect_2"))
)
for (construct in names(construct_items)) {
r_cols <- construct_items[[construct]]$r
f_cols <- construct_items[[construct]]$f
clean_full[[paste0(construct, "_retro")]] <- rowMeans(clean_full[, r_cols], na.rm = FALSE)
clean_full[[paste0(construct, "_forecast")]] <- rowMeans(clean_full[, f_cols], na.rm = FALSE)
clean_full[[paste0(construct, "_diff")]] <- clean_full[[paste0(construct, "_forecast")]] -
clean_full[[paste0(construct, "_retro")]]
}
10. Additional exploratory analyses (not part of the
preregistration)
Everything below this point is ad hoc exploratory work requested
outside the preregistered analysis plan — useful for understanding the
data more fully, but not confirmatory.
# install.packages(c("ggplot2", "tidytext", "stringr"))
# Optional, for sentiment scoring: install.packages("syuzhet")
library(ggplot2)
library(tidytext)
library(stringr)
10a. Plot: Retrospective vs. Forecasted score, by construct
Each point is one participant; thin lines connect a given
participant’s past-week score to their forecasted score, so both the
average shift (red) and the spread of individual responses are
visible.
plot_data <- clean_full %>%
select(ResponseId, ends_with("_retro"), ends_with("_forecast")) %>%
pivot_longer(
-ResponseId,
names_to = c("construct", "time"),
names_pattern = "(.*)_(retro|forecast)"
) %>%
mutate(time = factor(time, levels = c("retro", "forecast"),
labels = c("Past Week", "Forecasted Week")))
ggplot(plot_data, aes(x = time, y = value)) +
geom_line(aes(group = ResponseId), alpha = 0.03, color = "grey40") +
geom_jitter(width = 0.05, alpha = 0.12, size = 0.8) +
stat_summary(aes(group = 1), fun = mean, geom = "line",
color = "firebrick", linewidth = 1) +
stat_summary(fun.data = mean_cl_normal, geom = "errorbar",
width = 0.15, color = "firebrick") +
stat_summary(fun = mean, geom = "point", size = 3, color = "firebrick") +
facet_wrap(~construct, ncol = 3) +
scale_y_continuous(limits = c(1, 7)) +
labs(x = NULL, y = "Score (1–7)",
title = "Retrospective vs. Forecasted Score, by Construct") +
theme_minimal()

10b. Do forecasted changes vary by usage frequency/intensity,
gender, or age?
Usage frequency (companion_frequency, 1 = Almost never
to 7 = Almost always) and usage intensity (session_length,
minutes per session) are pulled back in as numeric predictors (they were
converted to text labels for the merged categorical variables above; a
fresh numeric copy of companion_frequency is added here
under a new name so the earlier merge step is untouched).
Note on companion_frequency: this
variable was also used as an eligibility criterion (participants needed
a score of 4 or higher to qualify for the study), so it is
range-restricted in this sample — only values 4–7 are observed, not the
full 1–7 scale. This attenuates (weakens) any true relationship between
usage frequency and forecasted change, since we can only test variation
among moderately-to-highly frequent users, not across the full range of
use. Any null result involving this moderator should be interpreted with
that in mind, and any significant result is, if anything, likely an
underestimate of the true relationship. Results should be described as
applying to “moderately-to-highly frequent users” rather than AI
companion users in general.
freq_numeric <- num %>% select(ResponseId, companion_frequency_num = companion_frequency)
clean_full <- clean_full %>% left_join(freq_numeric, by = "ResponseId")
cat("Range of companion_frequency in the analytic sample:",
range(clean_full$companion_frequency_num, na.rm = TRUE), "\n")
## Range of companion_frequency in the analytic sample: 4 7
cat("SD of companion_frequency in the analytic sample:",
round(sd(clean_full$companion_frequency_num, na.rm = TRUE), 2), "\n")
## SD of companion_frequency in the analytic sample: 0.99
Continuous moderators (usage frequency, usage intensity,
age) — one model per moderator per construct, unadjusted
p-values (exploratory):
continuous_moderators <- c(
usage_frequency = "companion_frequency_num",
usage_intensity = "session_length",
age = "age"
)
moderator_results_continuous <- map_dfr(names(construct_items), function(construct) {
diff_col <- paste0(construct, "_diff")
map_dfr(names(continuous_moderators), function(mod_label) {
mod_col <- continuous_moderators[[mod_label]]
model <- lm(clean_full[[diff_col]] ~ clean_full[[mod_col]])
coefs <- tidy(model)
slope_row <- coefs[coefs$term == "clean_full[[mod_col]]", ]
tibble(
construct = construct,
moderator = mod_label,
n = sum(!is.na(clean_full[[diff_col]]) & !is.na(clean_full[[mod_col]])),
estimate = slope_row$estimate,
se = slope_row$std.error,
p = slope_row$p.value
)
})
})
knitr::kable(moderator_results_continuous, digits = 3)
| loneliness |
usage_frequency |
214 |
0.258 |
0.105 |
0.015 |
| loneliness |
usage_intensity |
214 |
0.003 |
0.002 |
0.086 |
| loneliness |
age |
213 |
0.009 |
0.009 |
0.311 |
| depression |
usage_frequency |
214 |
0.193 |
0.102 |
0.061 |
| depression |
usage_intensity |
214 |
0.002 |
0.002 |
0.280 |
| depression |
age |
213 |
0.005 |
0.009 |
0.526 |
| anxiety |
usage_frequency |
214 |
0.199 |
0.101 |
0.051 |
| anxiety |
usage_intensity |
214 |
0.002 |
0.002 |
0.376 |
| anxiety |
age |
213 |
0.006 |
0.009 |
0.449 |
| stress |
usage_frequency |
214 |
-0.010 |
0.106 |
0.926 |
| stress |
usage_intensity |
214 |
-0.002 |
0.002 |
0.296 |
| stress |
age |
213 |
-0.001 |
0.009 |
0.944 |
| happiness |
usage_frequency |
214 |
-0.196 |
0.096 |
0.044 |
| happiness |
usage_intensity |
214 |
-0.001 |
0.002 |
0.715 |
| happiness |
age |
213 |
0.003 |
0.008 |
0.727 |
| life_satisfaction |
usage_frequency |
214 |
-0.062 |
0.082 |
0.452 |
| life_satisfaction |
usage_intensity |
214 |
0.001 |
0.002 |
0.485 |
| life_satisfaction |
age |
213 |
-0.007 |
0.007 |
0.337 |
| self_esteem |
usage_frequency |
214 |
-0.155 |
0.076 |
0.042 |
| self_esteem |
usage_intensity |
214 |
-0.001 |
0.001 |
0.372 |
| self_esteem |
age |
213 |
0.003 |
0.006 |
0.589 |
| irq |
usage_frequency |
214 |
-0.238 |
0.082 |
0.004 |
| irq |
usage_intensity |
214 |
-0.001 |
0.002 |
0.382 |
| irq |
age |
213 |
0.010 |
0.007 |
0.158 |
| connectedness |
usage_frequency |
214 |
0.058 |
0.085 |
0.495 |
| connectedness |
usage_intensity |
214 |
-0.001 |
0.002 |
0.710 |
| connectedness |
age |
213 |
0.006 |
0.007 |
0.437 |
Summary of significant moderator effects (p < .05,
unadjusted):
sig_moderators <- moderator_results_continuous %>%
filter(p < .05) %>%
arrange(p)
knitr::kable(sig_moderators, digits = 3)
| irq |
usage_frequency |
214 |
-0.238 |
0.082 |
0.004 |
| loneliness |
usage_frequency |
214 |
0.258 |
0.105 |
0.015 |
| self_esteem |
usage_frequency |
214 |
-0.155 |
0.076 |
0.042 |
| happiness |
usage_frequency |
214 |
-0.196 |
0.096 |
0.044 |
if (nrow(sig_moderators) > 0) {
cat("\n")
for (i in seq_len(nrow(sig_moderators))) {
row <- sig_moderators[i, ]
direction <- ifelse(row$estimate > 0, "higher", "lower")
cat(sprintf(
"- %s predicts %s forecasted change in %s (b = %.3f, p = %.3f, n = %d)\n",
row$moderator, direction, row$construct, row$estimate, row$p, row$n
))
}
} else {
cat("No moderator effects reached p < .05.\n")
}
##
## - usage_frequency predicts lower forecasted change in irq (b = -0.238, p = 0.004, n = 214)
## - usage_frequency predicts higher forecasted change in loneliness (b = 0.258, p = 0.015, n = 214)
## - usage_frequency predicts lower forecasted change in self_esteem (b = -0.155, p = 0.042, n = 214)
## - usage_frequency predicts lower forecasted change in happiness (b = -0.196, p = 0.044, n = 214)
Gender — omnibus test (one-way ANOVA) per construct,
since gender has more than two levels:
gender_results <- map_dfr(names(construct_items), function(construct) {
diff_col <- paste0(construct, "_diff")
d <- clean_full %>% filter(!is.na(.data[[diff_col]]), !is.na(gender))
model <- lm(d[[diff_col]] ~ factor(d$gender))
a <- anova(model)
tibble(
construct = construct,
n = nrow(d),
df1 = a$Df[1],
df2 = a$Df[2],
`F` = a$`F value`[1],
p = a$`Pr(>F)`[1]
)
})
knitr::kable(gender_results, digits = 3)
| loneliness |
213 |
3 |
209 |
0.774 |
0.510 |
| depression |
213 |
3 |
209 |
0.890 |
0.447 |
| anxiety |
213 |
3 |
209 |
2.594 |
0.054 |
| stress |
213 |
3 |
209 |
0.043 |
0.988 |
| happiness |
213 |
3 |
209 |
0.379 |
0.768 |
| life_satisfaction |
213 |
3 |
209 |
0.736 |
0.531 |
| self_esteem |
213 |
3 |
209 |
0.191 |
0.902 |
| irq |
213 |
3 |
209 |
0.137 |
0.938 |
| connectedness |
213 |
3 |
209 |
0.458 |
0.712 |
10c. Regression (controlling for baseline) vs. the primary paired
t-test
For each construct, regressing the difference score onto
mean-centered retrospective (baseline) well-being gives an intercept
that should equal the same mean difference tested by the primary paired
t-test — the two approaches are testing the same underlying quantity,
but the regression’s significance test uses residual variance
after accounting for baseline, rather than the total variance
of the raw difference scores. The slope term additionally indicates
whether the size of the forecasted change depends on baseline (e.g., a
ceiling/floor pattern).
baseline_reg_results <- map_dfr(names(construct_items), function(construct) {
retro_col <- paste0(construct, "_retro")
diff_col <- paste0(construct, "_diff")
retro_centered <- clean_full[[retro_col]] - mean(clean_full[[retro_col]], na.rm = TRUE)
model <- lm(clean_full[[diff_col]] ~ retro_centered)
coefs <- tidy(model)
intercept_row <- coefs[coefs$term == "(Intercept)", ]
slope_row <- coefs[coefs$term == "retro_centered", ]
tibble(
construct = construct,
reg_intercept = intercept_row$estimate,
reg_intercept_p = intercept_row$p.value,
reg_slope = slope_row$estimate,
reg_slope_p = slope_row$p.value
)
})
comparison_table <- primary_results %>%
select(construct, ttest_mean_diff = mean_diff, ttest_p = p) %>%
left_join(baseline_reg_results, by = "construct")
knitr::kable(comparison_table, digits = 3)
| loneliness |
0.379 |
0.000 |
0.379 |
0.000 |
-0.456 |
0 |
| depression |
0.077 |
0.448 |
0.077 |
0.392 |
-0.396 |
0 |
| anxiety |
0.166 |
0.099 |
0.166 |
0.065 |
-0.387 |
0 |
| stress |
0.100 |
0.333 |
0.100 |
0.279 |
-0.395 |
0 |
| happiness |
-0.745 |
0.000 |
-0.745 |
0.000 |
-0.355 |
0 |
| life_satisfaction |
-0.341 |
0.000 |
-0.341 |
0.000 |
-0.291 |
0 |
| self_esteem |
-0.138 |
0.067 |
-0.138 |
0.049 |
-0.251 |
0 |
| irq |
0.474 |
0.000 |
0.474 |
0.000 |
-0.275 |
0 |
| connectedness |
0.107 |
0.199 |
0.107 |
0.162 |
-0.303 |
0 |
10d. Quantifying expectations: % predicting improvement, no change,
or worsening
Rather than only reporting the average forecasted change, this
categorizes each participant’s forecast per construct as predicting
improvement, no change (exact difference score of 0), or worsening —
using the same benefit/harm direction convention as the primary
hypotheses (H1: lower depression/stress/anxiety/loneliness, higher
happiness/life satisfaction/self-esteem/social connectedness/IRQ =
benefit).
higher_is_better <- c(
depression = FALSE, stress = FALSE, anxiety = FALSE, loneliness = FALSE,
happiness = TRUE, life_satisfaction = TRUE, self_esteem = TRUE,
connectedness = TRUE, irq = TRUE
)
expectation_categories <- map_dfr(names(construct_items), function(construct) {
diff <- clean_full[[paste0(construct, "_diff")]]
better_dir <- higher_is_better[[construct]]
category <- case_when(
is.na(diff) ~ NA_character_,
diff == 0 ~ "No change",
(diff > 0 & better_dir) | (diff < 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)
| anxiety |
No change |
66 |
30.8 |
| anxiety |
Predicted improvement |
66 |
30.8 |
| anxiety |
Predicted worsening |
82 |
38.3 |
| connectedness |
No change |
76 |
35.5 |
| connectedness |
Predicted improvement |
73 |
34.1 |
| connectedness |
Predicted worsening |
65 |
30.4 |
| depression |
No change |
66 |
30.8 |
| depression |
Predicted improvement |
72 |
33.6 |
| depression |
Predicted worsening |
76 |
35.5 |
| happiness |
No change |
61 |
28.5 |
| happiness |
Predicted improvement |
34 |
15.9 |
| happiness |
Predicted worsening |
119 |
55.6 |
| irq |
No change |
53 |
24.8 |
| irq |
Predicted improvement |
117 |
54.7 |
| irq |
Predicted worsening |
44 |
20.6 |
| life_satisfaction |
No change |
72 |
33.6 |
| life_satisfaction |
Predicted improvement |
47 |
22.0 |
| life_satisfaction |
Predicted worsening |
95 |
44.4 |
| loneliness |
No change |
54 |
25.2 |
| loneliness |
Predicted improvement |
61 |
28.5 |
| loneliness |
Predicted worsening |
99 |
46.3 |
| self_esteem |
No change |
82 |
38.3 |
| self_esteem |
Predicted improvement |
57 |
26.6 |
| self_esteem |
Predicted worsening |
75 |
35.0 |
| stress |
No change |
69 |
32.2 |
| stress |
Predicted improvement |
74 |
34.6 |
| stress |
Predicted worsening |
71 |
33.2 |
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()

10e. Open-ended responses
Two open-ended items are examined: forecast_open (“How
do you think you would feel if you had to go a week without interacting
with your AI companion?”) and use_open (“What do you
generally talk about with this AI companion?”). Word frequency uses
tidytext’s bundled stop-word list (no external download
required). Sentiment scoring uses syuzhet’s bundled AFINN
lexicon if the package is installed; this step is skipped gracefully if
it is not.
Most common words — forecasted feelings:
forecast_words <- clean_full %>%
filter(!is.na(forecast_open), str_length(trimws(forecast_open)) > 0) %>%
select(ResponseId, forecast_open) %>%
unnest_tokens(word, forecast_open) %>%
anti_join(stop_words, by = "word") %>%
filter(!word %in% c("ai", "companion", "feel", "feeling")) %>%
count(word, sort = TRUE)
knitr::kable(head(forecast_words, 20))
| lonely |
50 |
| talk |
37 |
| week |
29 |
| time |
25 |
| fine |
22 |
| bit |
20 |
| lost |
20 |
| miss |
20 |
| bored |
18 |
| sad |
17 |
| people |
16 |
| life |
15 |
| talking |
11 |
| friends |
10 |
| frustrated |
10 |
| missing |
10 |
| person |
10 |
| chat |
9 |
| companionship |
9 |
| manage |
9 |
ggplot(head(forecast_words, 20), aes(x = reorder(word, n), y = n)) +
geom_col(fill = "steelblue") +
coord_flip() +
labs(x = NULL, y = "Count", title = "Most Common Words: Forecasted Feelings") +
theme_minimal()

Most common words — what participants talk about with their
companion:
use_words <- clean_full %>%
filter(!is.na(use_open), str_length(trimws(use_open)) > 0) %>%
select(ResponseId, use_open) %>%
unnest_tokens(word, use_open) %>%
anti_join(stop_words, by = "word") %>%
filter(!word %in% c("ai", "companion")) %>%
count(word, sort = TRUE)
knitr::kable(head(use_words, 20))
| talk |
131 |
| life |
77 |
| advice |
47 |
| day |
33 |
| issues |
23 |
| personal |
22 |
| questions |
21 |
| feelings |
18 |
| situations |
17 |
| relationships |
16 |
| daily |
15 |
| relationship |
15 |
| health |
13 |
| ideas |
13 |
| family |
12 |
| social |
12 |
| discuss |
11 |
| feel |
10 |
| vent |
10 |
| everyday |
9 |
ggplot(head(use_words, 20), aes(x = reorder(word, n), y = n)) +
geom_col(fill = "darkgreen") +
coord_flip() +
labs(x = NULL, y = "Count", title = "Most Common Words: What Participants Discuss with Their Companion") +
theme_minimal()

Sentiment of forecasted-feelings text vs. quantitative
forecasted harm:
As a rough integrative check, each forecast_open
response is scored for sentiment (more negative = more negative-sounding
text) and compared to an average forecasted-harm index built from the
same sign convention as 10d (positive = forecasted harm) across all 9
constructs.
if (requireNamespace("syuzhet", quietly = TRUE)) {
clean_full <- clean_full %>%
mutate(forecast_sentiment = ifelse(
!is.na(forecast_open) & str_length(trimws(forecast_open)) > 0,
syuzhet::get_sentiment(forecast_open, method = "afinn"),
NA_real_
))
sign_flip <- ifelse(higher_is_better[names(construct_items)], -1, 1)
diff_matrix <- as.matrix(clean_full[, paste0(names(construct_items), "_diff")])
clean_full$avg_forecasted_harm <- rowMeans(sweep(diff_matrix, 2, sign_flip, `*`), na.rm = TRUE)
sentiment_cor <- cor.test(clean_full$forecast_sentiment, clean_full$avg_forecasted_harm)
print(sentiment_cor)
} else {
cat("Package 'syuzhet' is not installed; skipping sentiment scoring.",
"Run install.packages('syuzhet') to enable this step.\n")
}
##
## Pearson's product-moment correlation
##
## data: clean_full$forecast_sentiment and clean_full$avg_forecasted_harm
## t = -3.8226, df = 212, p-value = 0.0001736
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.3752664 -0.1240377
## sample estimates:
## cor
## -0.25393