Abstract

Researchers working with heterogeneous experimental data frequently divide participants into subgroups. But this is typically done by hand, reintroducing the interpretive biases that experimental methods are meant to avoid. This talk advocates for a data-driven alternative: K-means clustering applied to participant response profiles.

I illustrate the approach with a truth-value judgment experiment on presupposition projection from quantifiers, which initially appeared to have failed: aggregate results did not conform to the (uncontroversial) predictions of any theory. Clustering revealed two subgroups, one using a non-linguistic response strategy and one patterning as predicted by the theoretical accounts. Bayesian distributional models fit per cluster confirm this and additionally support the variance hypothesis—that presupposition failure induces increased response variability due to Strawsonian squeamishness. Data-driven clustering, I argue, can rescue datasets that would otherwise be dismissed.

The experiment that I will present later is joint work, conceived and carried out together with Andrea Matticchio (UMass Amherst). All errors to do with the clustering are mine.

My goal for today will be a modest one: I want to show that it is possible to divide participants into groups, say because you find a bimodal distribution of responses, without having to make any of those categorization decisions yourself. That is, I will advocate for a process of determining participant clusters in a data-driven way, rather than in a researcher-intuition one.

Code

library(here)
library(tidyverse)
library(cluster)
library(ggalt)
library(ggrepel)
library(patchwork)
library(gghalves)
library(glue)
library(brms)
library(emmeans)
library(tidybayes)
library(bayestestR)
library(ggdist)
library(qrcode)

options(width = 10000)

source(here("scripts", "theme.R"))

Here are two participants from the experiment that I will be presenting later. Each of these is a representative of a different cluster. Disregarding for the moment what the labels on the axes mean, one cluster treats the undefined scenario just like the false and critical one, while in the other undefined appears to be its own category, being the only scenario that does not show an effect of negation.

Code

d <- read_csv(here("data", "all-projection.csv"), show_col_types = FALSE) %>%
    filter(trigger != "none") %>%
    select(id, scenario, negation, trigger, judgment) %>%
    mutate(
        id = as.numeric(as.factor(id)),
        negation = if_else(negation == "Alle", "without", "with"),
        judgment = judgment - 50,
        judgment = judgment / 25,
        scenario = gsub("undef", "undefined", scenario),
        cond = paste(scenario, negation, trigger, sep = "_")
    )

d %>%
    filter(id %in% c(10, 35)) %>%
    mutate(
        scenario = gsub("undefined", "undef.", scenario),
        scenario = gsub("critical", "crit.", scenario),
        scenario = fct_relevel(scenario, "true", "false", "undef."),
        negation = fct_relevel(negation, "without", "with"),
        id = paste("Participant", id)
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        color = negation,
        pch = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(fun = mean, geom = "point", size = 4) +
    guides(group = "none") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    coord_cartesian(clip = "off") +
    facet_wrap(~id) +
    theme(axis.text.x = element_markdown(family = "Cascadia Code"))

Over the course of this talk, I will show that it is possible to perform group assignments without having to involve the biases of the researcher, who is likely to be pursuing some kind of hypothesis herself.

The empirical domain we will be having a look at is presupposition projection, i.e., the question of how we are to predict the presuppositions of a complex sentence based on its parts.

We will start with a bit of background on presuppositions, though we do not need very much before we can dive into the experiment since we will leave out the theoretical hypothesis today.

Background: Presuppositions and truth-value intuitions

Natural language has expressions that trigger so-called presuppositions. These are prerequisites that must hold in order for the sentence to be meaningful to begin with. Most importantly for us, this amounts to truth-evaluability: an utterance that does not have its presuppositions met cannot be assigned a truth value (like true or false).

Compare believe and know below. The two sentences appear to differ with respect to whether the sentence conveys that Aditi actually owns a pet cockroach. The believe sentence comes with no such commitment, the know sentence however does.

1. Skye believes that Aditi owns a pet cockroach.
  \(\not\rightsquigarrow\) Aditi owns a pet cockroach.
2. Skye knows that Aditi owns a pet cockroach.
  \(\rightsquigarrow\) Aditi owns a pet cockroach.

If we explicitly negate this commitment, only the sentence with believe remains coherent. The minimal variant with know suffers from what is typically called presupposition failure and is felt to be degraded. Speakers very often have visceral reaction to presupposition failure, from hesitancy all the way to disgust. Strawson (1964) calls this reaction squeamishness.

1. Aditi does not own a pet cockroach but Skye believes that Aditi owns a pet cockroach.
2. #Aditi does not own a pet cockroach but Skye knows that Aditi owns a pet cockroach.

So know seems to be a presupposition trigger. There are other presupposition triggers, though:

1. Skye stopped smoking.
2. #Skye never smoked and they stopped smoking.
1. Skye smoked again.
2. #Skye never smoked and they smoked again.

Negation generally has no effect on the presuppositional requirement of an utterance. That is the negated and unnegated version of the same sentence involving a presupposition trigger have the same presuppositional requirement.

1. Skye stopped smoking.
  \(\rightsquigarrow\) Skye used to smoke.
2. Skye did not stop smoking.
  \(\rightsquigarrow\) Skye used to smoke.

By extension, a sentence suffering from presupposition failure without negation will generally still suffer from it when negated:

1. #Skye never smoked and they stopped smoking.
2. #Skye never smoked and they didn’t stop smoking.

For us, the main empirical data point will be about this particular property related to negation.

Projection from quantifiers over individuals: A failed experiment?

We, my collaborator Andrea Matticchio (UMass Amherst) and I, made use of a slider-based truth-value judgment task with the three labels completely false, neither nor, and completely true. This is a method that should closely mimic the well-known truth-value judgment task, allow for judgments indicating intuitions of presupposition failure via the “neither nor” option (Križ and Chemla 2015), and allow for the measurement of variability in judgments via the standard deviations associated with the condition means (see Sun et al. 2025 for discussion).

Code

ggplot() +
  geom_segment(aes(x = 0, xend = 1, y = 0, yend = 0), linewidth = 1) +
  geom_segment(aes(x = 0, xend = 0, y = -0.01, yend = 0.01), linewidth = 1) +
  geom_segment(aes(x = 1, xend = 1, y = -0.01, yend = 0.01), linewidth = 1) +
  geom_point(aes(x = 0.5, y = -0.01), shape = 24, size = 4, fill = "black") +
  annotate("text", x = 0, y = 0.05, label = "(completely false)", hjust = .5) +
  annotate("text", x = 0.5, y = 0.05,
           label = "(neither completely true nor completely false)",
           hjust = 0.5) +
  annotate("text", x = 1, y = 0.05, label = "(completely true)", hjust = .5) +
  xlim(-.1, 1.1) +
  ylim(-0.1, 0.1) +
  theme_void()

We will look at presupposition projection from quantifiers over individuals like every. We will have to set aside the important discussion of previous theoretical and experimental work on presupposition projection, but the interested reader is referred to Thalmann and Matticchio (2024);Thalmann and Matticchio (submitted) and the references therein.

Every woman canoed again.
Maybe presupposes: Every woman canoed before.

Here, we will mostly focus on the control conditions, one of which featured uncontroversial presupposition failure, no matter which theoretical framework is adopted.

The experiment was set up as follows:

\(2 \times 4\) design (48 participants, 48 items)
1. negation: without vs. with matrix negation (within items and participants)
2. scenario: true vs. false vs. undefined vs. critical (within items and participants)
3. Presupposition trigger as a pseudo-factor: again vs. stop (between items, 24 each)

Below is a sample item translated from the original German with the hard trigger again, plus the verbal context that instantiates the scenario manipulation. The visual context was informationally redundant and shown alongside the verbal context. The original item is presented in Section 6.1.

All women canoed again.
All women did not canoe again.
1. true:
  Maria canoed in the past and recently.
  Johanna canoed in the past and recently.
  Susanne canoed in the past and recently.
2. false:
  Maria canoed in the past and recently.
  Johanna canoed in the past and recently.
  Susanne canoed in the past but not recently.
3. undefined:
  Maria did not canoe in the past but did recently.
  Johanna did not canoe in the past and did not recently.
  Susanne canoed in the past and recently.
4. critical:
  Maria did not canoe in the past but did recently.
  Johanna canoed in the past and recently.
  Susanne canoed in the past but not recently.

The inclusion of the pseudo-factor was motivated by examples such as these. Arguably, these show that some presupposition triggers are less sensitive to presupposition failure contexts than others (Abusch 2002; Abusch 2010; Abrusán 2016). Though opinions on all of this are various, we will say here that this is due to local accommodation, which turns presupposition failure into mere falsity and which only applies to certain triggers like discover and stop (but not to others like regret and again).

1. I don’t know if Skye got Aditi a present, but if Aditi discovers that Skye did, she’ll give them a hug.
2. #I don’t know if Skye got Aditi a present, but if Aditi is happy that Skye did, she’ll give them a hug.
1. What did Skye buy? — Aditi discovered that they bought a guitar.
2. What did Skye buy? — #Aditi is happy that they bought a guitar.

That is, in our experiment, we have one hard trigger (again), which should not allow for local accommodation, and a soft trigger (stop), which according to received wisdom and some experimental results should (Charlow 2009; Chemla 2009; Fox 2013; Chen et al. 2022).

As I said before, I will not go into the theoretical details behind the hypotheses we’re testing here; see Thalmann and Matticchio (2024);Thalmann and Matticchio (submitted) for a full discussion (for quantifiers over worlds, i.e., attitude predicates). The key point is that the two theories we are comparing here come apart in their predictions for the critical scenario. People in the universal projection camp expect that this scenario should elicit intuitions of presupposition failure (i.e., judgments near “neither nor”) regardless of negation. By contrast, the strong Kleene theorist expects that negation should flip the judgment from “completely false” to “completely true” in this scenario.

Code

preds <- tribble(
    ~"Scenario" , ~"Negation" , ~"Judgment" , ~"who"                            ,
    "true"      , "without"   ,           2 , "*again*<br>Strong Kleene"        ,
    "true"      , "without"   ,           2 , "*again*<br>Universal Projection" ,
    "true"      , "without"   ,           2 , "*stop*<br>Local Accommodation"   ,
    "true"      , "with"      ,          -2 , "*again*<br>Strong Kleene"        ,
    "true"      , "with"      ,          -2 , "*again*<br>Universal Projection" ,
    "true"      , "with"      ,          -2 , "*stop*<br>Local Accommodation"   ,
    "false"     , "without"   ,          -2 , "*again*<br>Strong Kleene"        ,
    "false"     , "without"   ,          -2 , "*again*<br>Universal Projection" ,
    "false"     , "without"   ,          -2 , "*stop*<br>Local Accommodation"   ,
    "false"     , "with"      ,           2 , "*again*<br>Strong Kleene"        ,
    "false"     , "with"      ,           2 , "*again*<br>Universal Projection" ,
    "false"     , "with"      ,           2 , "*stop*<br>Local Accommodation"   ,
    "undef."    , "without"   ,           0 , "*again*<br>Strong Kleene"        ,
    "undef."    , "without"   ,           0 , "*again*<br>Universal Projection" ,
    "undef."    , "without"   ,          -2 , "*stop*<br>Local Accommodation"   ,
    "undef."    , "with"      ,           0 , "*again*<br>Strong Kleene"        ,
    "undef."    , "with"      ,           0 , "*again*<br>Universal Projection" ,
    "undef."    , "with"      ,           2 , "*stop*<br>Local Accommodation"   ,
    "crit."     , "without"   ,          -2 , "*again*<br>Strong Kleene"        ,
    "crit."     , "with"      ,           2 , "*again*<br>Strong Kleene"        ,
    "crit."     , "without"   ,           0 , "*again*<br>Universal Projection" ,
    "crit."     , "without"   ,          -2 , "*stop*<br>Local Accommodation"   ,
    "crit."     , "with"      ,           0 , "*again*<br>Universal Projection" ,
    "crit."     , "with"      ,           2 , "*stop*<br>Local Accommodation"   ,
)

preds %>%
    mutate(Negation = fct_relevel(Negation, "without")) %>%
    ggplot(aes(
        x = fct_inorder(Scenario),
        y = Judgment,
        color = Negation,
        shape = Negation,
        group = Negation
    )) +
    geom_line() +
    geom_point(size = 3) +
    scale_y_continuous(
        breaks = c(-2, 0, 2),
        labels = c("Compl. false", "Neither nor", "Compl. true")
    ) +
    labs(
        x = "Scenario",
        y = "Judgment"
    ) +
    coord_cartesian(clip = "off") +
    facet_wrap(~who) +
    theme(
        axis.text.y = element_markdown(hjust = 1),
        axis.text.x = element_markdown(family = "Cascadia Code", size = 8)
    )

Here we have the actual results:

Code

d %>%
    mutate(
        scenario = gsub("undefined", "undef.", scenario),
        scenario = gsub("critical", "crit.", scenario),
        scenario = fct_relevel(scenario, "true", "false", "undef."),
        negation = fct_relevel(negation, "without", "with"),
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        color = negation,
        pch = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = FALSE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = FALSE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(
        fun.data = mean_cl_normal,
        width = .15,
        alpha = .7,
        geom = "errorbar"
    ) +
    stat_summary(fun = mean, geom = "point", size = 4) +
    guides(group = "none") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    facet_wrap(~trigger) +
    theme(axis.text.x = element_markdown(family = "Cascadia Code"))

At first sight, this is just a complete mess: the two control conditions associated with the undefined scenario did not at all pattern as expected. We find neither scale-medial judgments nor an absence of an effect of negation.

In fact, there is an argument to be made that false, critical and undefined patterned exactly alike. Since this is not predicted by any of the theories on the market, this suggests that the experiment was a failure.

One plausible reason is that (some) participants at some point got tired of reading the linguistic material and simply interpreted the images. And since images have no presuppositions, no presupposition failure was detected.

But recall the plot from the outset which suggests that we do have participants who rated the items in accordance with our control expectations:

Code

d %>%
    filter(id %in% c(10, 35)) %>%
    mutate(
        scenario = gsub("undefined", "undef.", scenario),
        scenario = gsub("critical", "crit.", scenario),
        scenario = fct_relevel(scenario, "true", "false", "undef."),
        negation = fct_relevel(negation, "without", "with"),
        id = paste("Participant", id)
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        color = negation,
        pch = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(fun = mean, geom = "point", size = 4) +
    guides(group = "none") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    coord_cartesian(clip = "off") +
    facet_wrap(~id) +
    theme(axis.text.x = element_markdown(family = "Cascadia Code"))

I will now to offer a programmatic solution to assign participants to clusters. We will then plot our results again and do some statistical modelling to see whether the hypothesis can still be evaluated even though the experiment looks like a failure at first blush.

Clustering against the mess

When a dataset shows a bimodal distribution—that is, two distinct “humps” in the data rather than a single bell curve—it is often a signal that the sample contains two meaningfully different subgroups responding or behaving in fundamentally different ways. Though other types of response distributions can indicate multiple participant groups as well. Analysing the group as a whole in that situation can be misleading, as summary statistics like the mean may fall in the valley between the two humps, describing nobody in the sample particularly well.

Clustering is a natural response to this problem because it lets the data speak for itself: rather than imposing a single average profile on everyone, it allows the analysis to surface those underlying subgroups empirically. If a bimodal pattern reflects a genuine divide, clustering will tend to recover that structure and give each subgroup its own characterisation. This makes subsequent comparisons more meaningful, since differences observed between clusters are grounded in real variation in the data rather than being obscured by lumping everyone together.

One possible response to this is to simply group participants yourself. Here, I will pursue a different strategy, for the following reasons.

Assigning participants to groups by hand is tempting when patterns seem obvious, but it introduces a hard-to-avoid problem: the person doing the sorting already knows the data, and that knowledge inevitably shapes their decisions. This is a form of confirmation bias—groups end up reflecting what the researcher expected to find rather than what the data actually shows.

K-means clustering removes that subjectivity entirely by using a mathematical criterion to determine group membership, making the process transparent and reproducible. Anyone rerunning the same algorithm on the same data will get the same result, which is not something that can be said for hand-sorting.

Multidimensional scaling and K-means clustering

This will be the input for our clustering work, with participants as the rows and condition identifiers as the columns. Each non-identifier cell contains a condition mean, say the participants mean judgment for the critical scenario with negation and with the presupposition trigger again, in the column labelled critical_with_again.

Code

d_pivoted <- d %>%
    select(id, cond, judgment) %>%
    group_by(id, cond) %>%
    summarize(judgment = mean(judgment)) %>%
    pivot_wider(id_cols = id, names_from = cond, values_from = judgment)
d_pivoted

# A tibble: 48 × 17
# Groups:   id [48]
      id critical_with_again critical_with_stop critical_without_again critical_without_stop false_with_again false_with_stop false_without_again false_without_stop true_with_again true_with_stop true_without_again true_without_stop undefined_with_again undefined_with_stop undefined_without_again undefined_without_stop
   <dbl>               <dbl>              <dbl>                  <dbl>                 <dbl>            <dbl>           <dbl>               <dbl>              <dbl>           <dbl>          <dbl>              <dbl>             <dbl>                <dbl>               <dbl>                   <dbl>                  <dbl>
 1     1                1.8                1.76                 -1.92                  -2               1.71             0.6                -1.75              -1.64          -0.76          -1.96                1.67              1.59                0.227               1.37                    -1.61                  -1.83
 2     2                2                  2                    -2                     -2               2                2                  -2                 -2             -2             -2                   2                 2                   2                   2                       -2                     -2   
 3     3                1.75               1.59                 -1.89                  -0.48            1.77             1.83               -2                 -1.84          -1.35          -1.99                1.81              2                   1.61                1.82                    -1.73                  -1.41
 4     4                2                  2                    -1.97                  -2               2                2                  -1.97              -1.99          -2             -1.97                2                 2                   2                   0                       -1.96                  -1.97
 5     5                2                  2                    -2                     -2               2                2                  -2                 -2             -0.667         -2                   1.99              2                   0.667               2                       -2                     -2   
 6     6                2                  1.62                 -2                     -2               2                2                  -2                 -2             -1.2           -2                   2                 2                   0.28               NA                       -2                     -2   
 7     7                1.97               2                    -2                     -2               0.653            1.95               -2                 -1.96          -0.613         -0.693               0.64              1.99                0.667               0.613                   -1.97                  -2   
 8     8                2                  2                     0.667                 -2              -0.667            2                  -1.15              -2             -2             -2                   2                 2                   2                   2                       -2                     -2   
 9     9                2                  2                    -2                     -2               0.667            2                  -2                 -2             -0.667         -2                   2                 2                   1.52                0.667                   -2                     -0.92
10    10                2                  2                    -2                     -2               2                2                  -2                 -2             -2             -0.667               2                 2                   2                   2                       -2                     -2   
# ℹ 38 more rows

To identify natural groupings in the data, each subject is first represented as a point in a simplified two-dimensional space using a technique called Multidimensional Scaling (MDS). This step takes the full complexity of each subject’s measurements and “flattens” it into a map, where subjects with similar profiles end up close together and dissimilar subjects end up far apart.

So instead of having many condition means (a \(48 \times 17\) matrix as above), we reduce the complexity of the data to a \(48 \times 2\) matrix. Below, I only show a few of of the participants, since the scaled values are uninterpretable anyway.

Code

# Perform MDS analysis
n_cluster <- 3
mds_d <- cmdscale(dist(d_pivoted[, 2:length(d_pivoted)]))
head(mds_d, 10)

             [,1]         [,2]
 [1,]  0.43747216 -0.391105393
 [2,] -1.72309831 -0.400390137
 [3,] -0.78365437 -0.488091009
 [4,] -0.74803122  0.445352156
 [5,] -0.74281558 -0.188147979
 [6,]  0.13862523  0.657051192
 [7,]  0.44370312 -0.294300182
 [8,] -1.13920992 -1.257602271
 [9,] -0.18870645  0.034432487
[10,] -1.56538847 -0.708975981

Those points are then grouped using K-means clustering (Hartigan and Wong 1979), which partitions subjects into three clusters by iteratively finding the arrangement that minimises the distance between each point and the centre of its assigned group. The algorithm runs from many different starting positions to ensure a stable, reliable result. Each subject is ultimately assigned a cluster label, which is then joined back to the main dataset for further analysis.

Code

# Form clusters using K-means clustering
# changing the seed only affects the labelling of the clusters, not the sorting
set.seed(17)
kmeans_clusters <- kmeans(
    mds_d,
    iter.max = 5000,
    algorithm = "Hartigan-Wong",
    nstart = length(unique(d$id)),
    centers = n_cluster
)$cluster

mds_df <- as.data.frame(mds_d) %>%
    mutate(
        cluster = as.factor(kmeans_clusters),
        id = d_pivoted$id
    )
d_ac <- d %>%
    left_join(
        mds_df %>%
            select(id, cluster),
        by = "id"
    )

Here we have the clusters plotted with the the MDS axes. As you can see, we have three nicely separated clusters of participants. Though, of course, at this stage the clusters are not very informative because the two axes on the plot are not either.

Code

mds_df %>%
    ggplot(aes(
        x = V1,
        y = V2,
        color = cluster,
        fill = cluster,
        shape = cluster
    )) +
    geom_point(size = 1.5) +
    geom_encircle(
        alpha = .2,
        s_shape = 1,
        expand = 0,
        show.legend = FALSE,
        color = NA
    ) +
    geom_text_repel(
        data = filter(mds_df, cluster != "1"),
        aes(label = id),
        box.padding = 0.5,
        direction = "x",
        seed = 1234,
        size = 3,
        min.segment.length = 1,
        max.iter = 50000,
        max.time = 2,
        show.legend = FALSE
    ) +
    scale_shape_manual(values = c(16, 15, 17)) +
    labs(
        x = "Multidimensional scaling (MDS) for participant clusters; dimension 1",
        y = "MDS dimension 2",
        color = "Cluster",
        fill = "Cluster",
        shape = "Cluster"
    ) +
    scale_color_manual(values = colors5) +
    scale_fill_manual(values = colors5) +
    guides(shape = guide_legend(override.aes = list(size = 3))) +
    theme(legend.margin = margin(-10, 0, 0, 0, "pt"))

To remedy this, we can plot the descriptive results again, but separated by cluster.

Code

d_ac %>%
    arrange(cluster) %>%
    mutate(
        scenario = gsub("undefined", "undef.", scenario),
        scenario = gsub("critical", "crit.", scenario),
        scenario = fct_relevel(scenario, "true", "false", "undef."),
        negation = fct_relevel(negation, "without"),
        trigger = paste0("*", trigger, "*"),
        color = colors5[as.numeric(cluster)],
        cluster = glue("<i style='color:{color}'>**Cluster {cluster}**</i>"),
        cluster = fct_inorder(cluster)
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        pch = negation,
        color = negation,
        fill = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(fun.data = mean_cl_normal, geom = "errorbar", width = .1) +
    stat_summary(fun = mean, geom = "point", size = 2) +
    facet_wrap(trigger ~ cluster, ncol = n_cluster) +
    coord_cartesian(clip = "off") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        fill = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    guides(shape = guide_legend(override.aes = list(size = 4))) +
    theme(
        legend.margin = margin(-25, 0, 0, 0, "pt"),
        axis.text.x = element_markdown(family = "Cascadia Code", size = 7)
    )

Here, we now have two patterns that are easy to interpret:

Cluster 1: these participants used an interpretive strategy that mitigated presupposition failure, resulting in responses that do not align with theoretical predictions.¹
Cluster 2: these participants, on the other hand, have the characteristic feature of presupposition failure for the undefined scenario: scale-medial judgments and no (meaningful) effect of negation.

Since Cluster 3 consists of only two participants who simply misunderstood the task according to their judgments, we will leave this cluster aside from here on.

Clustered models

At this stage, we are finally in a position to do some statistical analysis that will allow us to evaluate the hypotheses I showed you earlier. To do this, we will first fit two models, one for each cluster. We will fit so-called distributional models. A distributional model like the one below estimates not just the average value of a response (the mean) but also how much variability surrounds that average, giving a fuller picture of how scores are spread across individuals or conditions.

1. Y ~ negation * scenario * trigger +
  (1 + scenario * negation * trigger | id) +
  (1 + scenario * negation | item)
2. \(\log \sigma\) ~ negation * scenario + trigger

We will return to the looming question why it makes sense to estimate not just means but also by-condition standard deviations for this type of experiment and empirical domain in a moment.

All the models here are Bayesian mixed models, fit with the brms package (Bürkner 2018; Bürkner 2021), that allow for the parallel estimation of the location and scale parameters of the judgment distributions. The priors are all weakly informative, biasing the estimation towards smaller effects centered on 0, but not ruling out larger ones.

Code

# run a bayesian model for each of the two big clusters (mixed distributional)
cluster2 <- d_ac %>%
    filter(cluster == "2") %>%
    pull(id) %>%
    unique()
cluster1 <- d_ac %>%
    filter(cluster == "1") %>%
    pull(id) %>%
    unique()
d <- read_csv(here("data", "all-projection.csv"), show_col_types = FALSE) %>%
    filter(trigger != "none") %>%
    mutate(
        id = as.numeric(as.factor(id)),
        negation = if_else(negation == "Alle", "without", "with"),
        judgment = judgment - 50,
        judgment = judgment / 25,
        scenario = fct_relevel(scenario, "false", "undef", "critical", "true"),
        negation = fct_relevel(negation, "with"),
        trigger = fct_relevel(trigger, "stop", "again")
    )

# function to have nicer names for sum coding
named.contr.sum <- function(x, ...) {
    if (is.factor(x)) {
        x <- levels(x)
    } else if (is.numeric(x) & length(x) == 1L) {
        stop("cannot create names with integer value. Pass factor levels")
    }
    x <- contr.sum(x, ...)
    colnames(x) <- paste0(rownames(x)[1:(nrow(x) - 1)], "_v_GM")
    return(x)
}

contrasts(d$negation) <- named.contr.sum(levels(d$negation))
contrasts(d$scenario) <- named.contr.sum(levels(d$scenario))
contrasts(d$trigger) <- named.contr.sum(levels(d$trigger))

mod_clust2 <- brm(
    bf(
        judgment ~ negation *
            scenario *
            trigger +
            (1 + scenario * negation * trigger | id) +
            (1 + scenario * negation | item),
        # note that these are modelled on the log scale; use exp(sigma) for interpretation
        sigma ~ negation * scenario + trigger
    ),
    data = d %>%
        filter(id %in% cluster2),
    family = gaussian(),
    prior = c(
        prior(normal(0, .5), class = Intercept),
        prior(normal(0, .5), class = Intercept, dpar = sigma),
        prior(normal(0, 1), class = b),
        prior(normal(0, .25), class = b, dpar = sigma),
        prior(normal(0, .5), class = sd),
        prior(lkj(2), class = cor)
    ),
    iter = 10000,
    seed = 1234,
    cores = 4,
    control = list(max_treedepth = 12, adapt_delta = .9),
    silent = 2,
    refresh = 0,
    file = here("models", "allproj-clust2"),
    chains = 4
)

mod_clust1 <- update(
    mod_clust2,
    newdata = d %>%
        filter(id %in% cluster1),
    file = here("models", "allproj-clust1"),
    cores = 4
)

Code

fixef(mod_clust1)

                                                               Estimate   Est.Error          Q2.5         Q97.5
Intercept                                                -0.03482270338 0.019468091 -0.0728327329  0.0032980864
sigma_Intercept                                          -0.48412459079 0.020677547 -0.5250287508 -0.4438489537
negationwith_v_GM                                         0.85798842836 0.028425072  0.8028044980  0.9142763455
scenariofalse_v_GM                                       -0.00853866718 0.030604888 -0.0682464662  0.0513549360
scenarioundef_v_GM                                       -0.15240036258 0.038780284 -0.2281891237 -0.0762278011
scenariocritical_v_GM                                     0.02326395306 0.028263179 -0.0320355760  0.0787300360
triggerstop_v_GM                                          0.00490369323 0.019091086 -0.0325049068  0.0424121049
negationwith_v_GM:scenariofalse_v_GM                      0.95944746199 0.039379865  0.8818812641  1.0363768371
negationwith_v_GM:scenarioundef_v_GM                      0.77045228407 0.043060362  0.6863742266  0.8552681901
negationwith_v_GM:scenariocritical_v_GM                   0.96487100472 0.032694029  0.9014655135  1.0297377991
negationwith_v_GM:triggerstop_v_GM                       -0.00512090274 0.026294995 -0.0568820434  0.0465689599
scenariofalse_v_GM:triggerstop_v_GM                       0.00589073302 0.029777318 -0.0527610289  0.0635717609
scenarioundef_v_GM:triggerstop_v_GM                       0.03126677468 0.039514617 -0.0467212965  0.1078831917
scenariocritical_v_GM:triggerstop_v_GM                   -0.00300858230 0.028887516 -0.0598671350  0.0533575294
negationwith_v_GM:scenariofalse_v_GM:triggerstop_v_GM     0.00062711479 0.038177870 -0.0744839316  0.0747781597
negationwith_v_GM:scenarioundef_v_GM:triggerstop_v_GM    -0.02403183654 0.041613822 -0.1071404102  0.0575671046
negationwith_v_GM:scenariocritical_v_GM:triggerstop_v_GM  0.07013321412 0.031724409  0.0072624919  0.1321725487
sigma_negationwith_v_GM                                   0.29307812966 0.019885062  0.2540446648  0.3320548440
sigma_scenariofalse_v_GM                                 -0.05442324849 0.035279857 -0.1233092266  0.0154760795
sigma_scenarioundef_v_GM                                  0.33628783949 0.034853990  0.2674632263  0.4040476521
sigma_scenariocritical_v_GM                              -0.11402031059 0.035891212 -0.1841620169 -0.0439397854
sigma_triggerstop_v_GM                                   -0.10146772527 0.020604653 -0.1415265566 -0.0607824972
sigma_negationwith_v_GM:scenariofalse_v_GM               -0.08580337268 0.035140386 -0.1543373759 -0.0170684267
sigma_negationwith_v_GM:scenarioundef_v_GM                0.00463533075 0.032206310 -0.0588627491  0.0675039011
sigma_negationwith_v_GM:scenariocritical_v_GM            -0.17995582245 0.035147020 -0.2477210390 -0.1106374708

Code

fixef(mod_clust2)

                                                              Estimate   Est.Error         Q2.5        Q97.5
Intercept                                                -0.2752118699 0.101625762 -0.475185916 -0.069490322
sigma_Intercept                                          -0.7727350104 0.057385445 -0.880294407 -0.653792401
negationwith_v_GM                                         0.4017962548 0.090488899  0.218380601  0.579781700
scenariofalse_v_GM                                        0.1167303265 0.084023581 -0.049439881  0.280717169
scenarioundef_v_GM                                       -0.4761238530 0.125838817 -0.724731204 -0.229224154
scenariocritical_v_GM                                     0.0458817326 0.085322894 -0.121368544  0.216403511
triggerstop_v_GM                                          0.0173641859 0.051827545 -0.084486336  0.120660356
negationwith_v_GM:scenariofalse_v_GM                      1.1745321353 0.124339477  0.920755253  1.418767183
negationwith_v_GM:scenarioundef_v_GM                      0.0624340054 0.118638705 -0.167432113  0.301159734
negationwith_v_GM:scenariocritical_v_GM                   1.0413541346 0.114476002  0.813246213  1.266281661
negationwith_v_GM:triggerstop_v_GM                        0.0098048781 0.050724913 -0.091203185  0.110114434
scenariofalse_v_GM:triggerstop_v_GM                       0.0096241746 0.077364508 -0.146220810  0.160343089
scenarioundef_v_GM:triggerstop_v_GM                       0.0263754458 0.095618154 -0.163910724  0.213761096
scenariocritical_v_GM:triggerstop_v_GM                    0.0870254640 0.074249516 -0.056641867  0.235976630
negationwith_v_GM:scenariofalse_v_GM:triggerstop_v_GM     0.0556093500 0.076671445 -0.095619680  0.205692088
negationwith_v_GM:scenarioundef_v_GM:triggerstop_v_GM    -0.0632020707 0.112010892 -0.282898800  0.160153018
negationwith_v_GM:scenariocritical_v_GM:triggerstop_v_GM  0.0477691601 0.072993573 -0.096623005  0.190009837
sigma_negationwith_v_GM                                   0.6675148881 0.055892626  0.551095628  0.771811555
sigma_scenariofalse_v_GM                                  0.4584358863 0.079874933  0.302900908  0.615104398
sigma_scenarioundef_v_GM                                  0.8859535406 0.079924173  0.727670611  1.038526235
sigma_scenariocritical_v_GM                               0.2040352258 0.082989079  0.037408195  0.365200001
sigma_triggerstop_v_GM                                   -0.1241888892 0.045892546 -0.214687769 -0.034106113
sigma_negationwith_v_GM:scenariofalse_v_GM               -0.3318531029 0.078776570 -0.484152339 -0.173817215
sigma_negationwith_v_GM:scenarioundef_v_GM               -0.5133761765 0.079915275 -0.666518832 -0.353074836
sigma_negationwith_v_GM:scenariocritical_v_GM            -0.2743313634 0.084376004 -0.437330428 -0.107981606

Granted, this is incredibly hard to interpret.

Below, we show the estimated distributions for the conditions in the experiment, based on the estimated marginal means associated with the model parameter (emmeans package, Lenth 2024). In the second cluster, the undefined scenario patterns exactly as predicted, with judgments near “neither nor” and no effect of negation.

Code

p_post <- function(mod, cluster) {
    spreads_crit <- mod %>%
        emmeans(specs = c("scenario", "negation", "trigger")) %>%
        gather_emmeans_draws() %>%
        ungroup() %>%
        mutate(
            cond = paste0(
                scenario,
                "-",
                negation
            ),
            .before = scenario
        )

    spreads_crit %>%
        mutate(
            trigger = paste0("*", trigger, "*"),
            trigger = fct_relevel(trigger, "*again*", "*stop*")
        ) %>%
        ggplot(aes(
            x = fct_reorder(cond, .value, .desc = TRUE),
            y = .value,
            color = colors[1],
            fill = colors[1]
        )) +
        geom_hline(
            yintercept = 0,
            color = colors[2],
            lty = "dashed",
            alpha = .5
        ) +
        stat_slab(
            aes(fill = after_stat(level)),
            position = "dodgejust",
            point_interval = mean_hdi,
            .width = c(.95, .999),
            limits = c(.15, .85),
            color = "white",
            normalize = "xy",
            scale = .8,
            n = 5000,
            linewidth = .1
        ) +
        stat_pointinterval(
            position = position_dodge(width = .4, preserve = "single"),
            .width = c(.5, .9, .95),
            limits = c(.15, .85),
            point_size = 1.5,
            fill = "white",
            pch = 22,
            stroke = .1,
            point_interval = mean_qi
        ) +
        guides(
            fill_ramp = "none",
            fill = "none",
            pch = "none",
            color = "none"
        ) +
        scale_x_discrete(expand = c(0, 0)) +
        scale_y_continuous(
            limits = c(-2.2, 2.2),
            breaks = c(-2, -1, 0, 1, 2),
            labels = c("-2\nC. false", "-1\n", "0\nNeither nor", "1\n", "2\nC. true")
        ) +
        scale_fill_manual(
            values = c(alpha(colors[2], .7), alpha(colors[1], .7))
        ) +
        coord_flip(clip = "off") +
        theme(
            panel.grid.major.y = element_blank(),
            panel.grid.minor.y = element_blank(),
            axis.text.y = element_markdown(hjust = 1),
            axis.title.x = element_markdown(),
            axis.text.x = element_markdown()
        ) +
        labs(
            x = "Parameter",
            subtitle = glue(
                "**<i style='color:{colors5[as.numeric(cluster)]}'> Cluster {cluster}**</i>"
            ),
            y = "Posterior Marginal Means \u00B1HDI<sub>95</sub> \u00B1ETI<sub>50;90;95</sub>"
        ) +
        facet_wrap(~trigger, scales = "fixed")
}

p_post(mod_clust1, "1") /
    p_post(mod_clust2, "2")

To give an intermediate summary, we have seen so far that despite the desperation-inducing descriptive plot of the experiment we still managed to arrive at a sensible interpretation of the data because of clustering. It turns out that there was only one (if dominating) subgroup of participants that somehow interpreted the images rather than the linguistic stimuli. The other—identified based on clustering—behaved exactly as predicted by all theories of presuppositions that predict a local presupposition.

In the next step, we will go a bit further and see whether the second cluster also features another diagnostic of presuppositions that goes beyond the missing effect of negation and instead targets the standard deviations we had our distributional model estimate.

A detour: The variance hypothesis

There is also a second pattern in the data that I want to go through briefly. But before we do so, recall from the background section that participants hate presupposition failure. I call this Strawsonian squeamishness.

Strawsonian squeamishness (based on Strawson 1964)
An utterance suffering from presupposition failure and thereby lacking a classical truth value will cause interlocutors to experience a sense of conversational unease or breakdown—manifested as hesitation, uncertainty, or reluctance to evaluate the utterance for truth.

In previous work, I argue that squeamishness in experimental participants is detectable via increased standard deviations for those conditions (Thalmann and Matticchio 2024; Thalmann in progress; Thalmann and Matticchio submitted); I call this the variance hypothesis.

Variance hypothesis (Thalmann in progress)
Experimental conditions involving presupposition failure will exhibit greater response variability, i.e., higher standard deviations, than conditions that are evaluable for truth, due to Strawsonian squeamishness.

Because the models that I fitted actually estimate standard deviation coefficients as well as means, we can now have a look at whether the variance hypothesis holds for our data as well. The plot below gives the posterior marginal standard deviations for each of the conditions in the two clusters.

Code

p_post_sigma <- function(mod, cluster) {
    spreads_sigma <- mod %>%
        emmeans(
            specs = c("scenario", "negation", "trigger"),
            dpar = "sigma"
        ) %>%
        gather_emmeans_draws() %>%
        ungroup() %>%
        mutate(
            cond = paste0(
                scenario,
                "-",
                negation
            ),
            .before = scenario
        )

    spreads_sigma %>%
        mutate(
            trigger = paste0("*", trigger, "*"),
            trigger = fct_relevel(trigger, "*again*", "*stop*")
        ) %>%
        ggplot(aes(
            x = fct_reorder(cond, .value, .desc = FALSE),
            y = exp(.value),
            color = colors[1],
            fill = colors[1]
        )) +
        stat_slab(
            aes(fill = after_stat(level)),
            position = "dodgejust",
            point_interval = mean_hdi,
            .width = c(.95, .999),
            limits = c(.15, .85),
            color = "white",
            scale = .8,
            normalize = "xy",
            n = 5000,
            linewidth = .1
        ) +
        stat_pointinterval(
            position = position_dodge(width = .4, preserve = "single"),
            .width = c(.5, .9, .95),
            limits = c(.15, .85),
            point_size = 1.5,
            fill = "white",
            pch = 22,
            stroke = .1,
            point_interval = mean_qi
        ) +
        guides(
            fill_ramp = "none",
            fill = "none",
            pch = "none",
            color = "none"
        ) +
        scale_x_discrete(expand = c(0, 0)) +
        scale_fill_manual(
            values = c(alpha(colors[2], .7), alpha(colors[1], .7))
        ) +
        scale_y_continuous(limits = c(0, NA)) +
        coord_flip(clip = "off") +
        theme(
            panel.grid.major.y = element_blank(),
            panel.grid.minor.y = element_blank(),
            axis.text.y = element_markdown(hjust = 1),
            axis.title.x = element_markdown(),
            axis.text.x = element_markdown()
        ) +
        labs(
            x = "Parameter",
            subtitle = glue(
                "**<i style='color:{colors5[as.numeric(cluster)]}'> Cluster {cluster}**</i>"
            ),
            y = "Posterior Marginal Standard Deviations \u00B1HDI<sub>95</sub> \u00B1ETI<sub>50;90;95</sub>"
        ) +
        facet_wrap(~trigger, scales = "fixed")
}
p_post_sigma(mod_clust1, "1") /
    p_post_sigma(mod_clust2, "2")

As we can see, the undefined scenario, which is the only one that should be associated with presupposition failure in terms of its location estimate, is also the one associated with the largest scale estimate. This is consistent with the variance hypothesis, which predicts that intuitions of presupposition failure should be associated with increased variability in judgments due to Strawsonian squeamishness.

Some actual hypothesis testing

To do our hypothesis testing in a Bayesian framework, we cannot rely on \(p\)-values. Instead, we will use Bayes Factors and ROPEs.

A Bayes Factor is a way of quantifying the evidence that the data provide for one hypothesis over another; for example, that two conditions genuinely differ versus that they do not. Rather than asking whether an effect is statistically significant in the traditional sense, it asks how much more probable the observed data are under one model compared to another. Working on a logarithmic scale is convenient because it turns this ratio into a difference, making large and small values easier to handle numerically, and because symmetric positive and negative values correspond naturally to evidence for and against an effect respectively.

A Region of Practical Equivalence (ROPE) addresses a limitation of asking simply whether an effect exists at all: in practice, a tiny difference between conditions may be statistically detectable but too small to matter. The ROPE defines a narrow band around zero—representing effect sizes so small they are considered practically meaningless—and the analysis then evaluates whether the evidence supports an effect that falls outside this band. Combining this with a Bayes Factor means the conclusion is not just about whether conditions differ, but whether they differ by enough to be worth caring about, which is a more useful question.

Code

p_bfs <- function(model, cluster) {
    mod_max_prior <- update(
        model,
        sample_prior = "only",
        cores = 4,
        silent = 2,
        refresh = 0,
        file = here("models", paste0("allproj_prior_clust", cluster))
    )
    # the parameters we want to test; i.e., without trigger distinction
    params <- list(
        scenario = c("true", "false", "undef", "critical"),
        negation = c("with", "without"),
        trigger = NA
    )

    # get the reference grids for the posterior and the prior
    mod_max_rg <- ref_grid(model, at = params)
    mod_max_prior_rg <- ref_grid(mod_max_prior, at = params)

    mod_max_ropes <- bayesfactor_rope(
        pairs(mod_max_rg),
        prior = pairs(mod_max_prior_rg),
        effects = "fixed"
    )

    p_bfs_rope <- mod_max_ropes %>%
        as_tibble() %>%
        mutate(
            Parameter = gsub(" NA", "", contrast),
            # reorder the words so that the order matches the other plot
            Parameter = paste(
                stringr::word(Parameter, 2),
                stringr::word(Parameter, 1),
                "vs.",
                stringr::word(Parameter, 5),
                stringr::word(Parameter, 4)
            ),
            interpretation = effectsize::interpret_bf(
                log_BF,
                log = TRUE,
                rules = "raftery1995"
            ),
            sub_zero = if_else(log_BF > 0, "no", "yes")
        ) %>%
        separate(
            Parameter,
            into = c("left", "right"),
            sep = " vs. ",
            remove = FALSE
        ) %>%
        filter(
            left == "undef with" & right == "undef without" |
            left == "true with" & right == "true without" |
            left == "false with" & right == "false without" |
            left == "critical with" & right == "critical without"
        ) %>%
        mutate(
            Parameter = gsub("> ", ">-", Parameter)
        ) %>%
        ggplot(aes(
            y = fct_reorder(Parameter, log_BF),
            x = log_BF,
            color = sub_zero,
            pch = sub_zero
        )) +
        geom_point(size = 3) +
        scale_x_continuous(limits = c(0, NA)) +
        guides(color = "none", pch = "none") +
        labs(
            subtitle = glue(
                "**<i style='color:{colors5[cluster]}'> Cluster {cluster}**</i>"
            ),
            y = "Parameter contrasts",
            x = "BF<sub>log</sub> for Region of Practical Equivalence (ROPE): 0 \u00B1 .1 * SD(Y)"
        ) +
        theme(
            axis.title.x = element_markdown(),
            axis.text.y = element_markdown(hjust = 1)
        )

    # standard deviation: get the reference grids for the posterior and the prior
    mod_max_rg_sigma <- ref_grid(model, dpar = "sigma", at = params)
    mod_max_prior_rg_sigma <- ref_grid(
        mod_max_prior,
        dpar = "sigma",
        at = params
    )

    # compute the bfs on the standard deviations
    mod_max_bfs_sigma <- bayesfactor_rope(
        pairs(mod_max_rg_sigma),
        prior = pairs(mod_max_prior_rg_sigma),
        effects = "fixed"
    )

    p_bfs_rope_sigma <- mod_max_bfs_sigma %>%
        as_tibble() %>%
        mutate(
            Parameter = gsub(" NA", "", contrast),
            # reorder the words so that the order matches the other plot
            Parameter = paste(
                stringr::word(Parameter, 2),
                stringr::word(Parameter, 1),
                "vs.",
                stringr::word(Parameter, 5),
                stringr::word(Parameter, 4)
            ),
            interpretation = effectsize::interpret_bf(
                log_BF,
                log = TRUE,
                rules = "raftery1995"
            ),
            sub_zero = if_else(log_BF > 0, "no", "yes")
        ) %>%
        separate(
            Parameter,
            into = c("left", "right"),
            sep = " vs. ",
            remove = FALSE
        ) %>%
        filter(
            left == "undef with" & right == "undef without" |
            left == "true with" & right == "true without" |
            left == "false with" & right == "false without" |
            left == "critical with" & right == "critical without"
        ) %>%
        mutate(
            Parameter = gsub("> ", ">-", Parameter)
        ) %>%
        ggplot(aes(
            y = fct_reorder(Parameter, log_BF),
            x = log_BF,
            color = sub_zero,
            pch = sub_zero
        )) +
        scale_x_continuous(limits = c(0, NA)) +
        geom_point(size = 3) +
        guides(color = "none", pch = "none") +
        labs(
            subtitle = glue(
                "**<i style='color:{colors5[cluster]}'> Cluster {cluster}**</i>"
            ),
            y = "Parameter contrasts, Sigmas",
            x = "BF<sub>log</sub> for Region of Practical Equivalence (ROPE): 0 \u00B1 .1 * SD(Y)"
        ) +
        theme(
            axis.title.x = element_markdown(),
            axis.text.y = element_markdown(hjust = 1)
        )

    list(p_bfs_rope, p_bfs_rope_sigma)
}
clust1_bfs <- p_bfs(mod_clust1, cluster = 1)
clust2_bfs <- p_bfs(mod_clust2, cluster = 2)

The plot below compares the condition means, i.e., the location estimates. For the second cluster, it indicates that the effect for negation does not have strong support in the undefined condition, as expected under an assumption of presupposition failure. With the other conditions, we do find a large effect of negation, suggesting (as predicted) that no presupposition failure occurred.

Code

clust1_bfs[[1]] /
    clust2_bfs[[1]] +
    plot_layout(axis_titles = "collect")

With the standard deviations, we again find that while the presence of negation increases standard deviations in the other conditions, there is really no evidence for the undefined scenario to assume so.

Code

clust1_bfs[[2]] /
    clust2_bfs[[2]] +
    plot_layout(axis_titles = "collect")

Though tentative (likely because of the small amount of data), our models post-clustering show that there does exist a subgroup of participants who did not fail the control conditions. In other words, they were sensitive to the manipulation of presupposition failure in the experiment.

This is a large improvement over the impression that the descriptive data pre-clustering gave us initially, where we felt forced to consider the experiment failed.

Conclusion

We have seen that in the face of a mixture of participant subgroups, K-means clustering is an easy and reasonable objective way to group participants. In our case, we had one group that had truly uninterpretable results, likely because of a nonlinguistic response strategy. The other group, instead, patterned as expected for the controls, allowing for further analysis targeting the initial hypothesis behind the experiment. All in all then, clustering restored the usefulness of a dataset that initially looked like it indicated a failed experiment.

Below, I go through a couple of caveats.

First, I never justified why I created \(3\) cluster as opposed to any other number. Since my point above was that clustering comes with some degree of objectivity, justifying that type of decision in an objective way too would have been necessary. For reasons of time, I moved that discussion and justification to the appendix. The interested reader can have a look at Section 6.2. In Section 6.4 I show what would have happened if I had only had two clusters (hint: nothing much different).

Second, it is possible to use the clusters differently than I did here. In particular, it might be worthwhile to not compute separate models for each cluster but to instead add cluster membership as a predictor in a single model. This way, one could estimate all of the effects in a single model, which would improve the estimation of by-item effects, since more data would be available.

There are also model improvements to be made. In particular, the gaussian family that I used is probably inappropriate. Two candidates are the skew_normal family (which I explored in my talk at the DGfS this year and which estimates skew parameters alongside means and standard deviations) and the Beta family. Both of these allow to take into account the boundedness of the scale in some way, which the gaussian family simply cannot. I leave this for later because fitting these models takes a lot of time and the point today was mainly about the clustering algorithm.

One question that one always has to answer when it comes to clustering is what data the clustering algorithm is to have access to. Here, I opted for classical condition means, but there are good arguments to consider location estimates that are better estimates than means, i.e., those that are estimated with shrinkage. For an explanation of what shrinkage is and why shrinkage makes means better for clustering participants, I refer the reader to a recent blog post by Nathaniel Haines.

Lastly, there is an issue with the variance hypothesis portion of the talk. As I show on my blog, standard deviation estimates are biased with a small number of observations, leading to a systematic underestimation of standard deviation coefficients. Here, because we only analyzed a subset of the experimental participants for each model, it is likely that our estimates for the standard deviations are too small in some conditions. This might explain why the false and undefined scenarios with negation showed comparable variance estimates in the second cluster. This is not a pattern that was found in Thalmann (in progress);Thalmann and Matticchio (submitted);Thalmann and Matticchio (2024), where the false showed lower standard deviation and where considerably more data were available. That experiment and the descriptive results are shown in Section 6.5.

You can find this handout here:

References

Abrusán, Márta. 2016. Presupposition cancellation: Explaining the “soft–hard” trigger distinction. Natural Language Semantics 24. Springer: 165–202. https://doi.org/10.1007/s11050-016-9122-7.

Abusch, Dorit. 2002. Lexical alternatives as a source of pragmatic presuppositions. In Proceedings of Semantics and Linguistic Theory (SALT) 12, ed. Brendan Jackson, 1–19. https://doi.org/10.3765/salt.v12i0.2867.

Abusch, Dorit. 2010. Presupposition triggering from alternatives. Journal of Semantics 27. Oxford University Press: 37–80. https://doi.org/10.1093/jos/ffp009.

Bürkner, Paul-Christian. 2018. Advanced bayesian multilevel modeling with the R package brms. The R Journal 10: 395–411. https://doi.org/10.32614/RJ-2018-017.

Bürkner, Paul-Christian. 2021. Bayesian item response modeling in R with brms and stan. Journal of Statistical Software 100: 1–54. https://doi.org/10.18637/jss.v100.i05.

Charlow, Simon. 2009. “Strong” predicative presuppositional objects. In Proceedings of ESSLLI 2009.

Chemla, Emmanuel. 2009. Presuppositions of quantified sentences: Experimental data. Natural Language Semantics 17. Springer: 299–340. https://doi.org/10.1007/s11050-009-9043-9.

Chen, Yuqiu, Maik Thalmann, and Mailin Antomo. 2022. Presupposition triggers and (not-)at-issueness. Insights from language acquisition into the soft-hard distinction. Journal of Pragmatics 199: 21–46. https://doi.org/10.1016/j.pragma.2022.06.014.

Fox, Danny. 2013. Presupposition projection from quantificational sentences: Trivalence, local accommodation, and presupposition strengthening. In From grammar to meaning, ed. Ivano Caponigro and Carlo Cecchetto, 201–232. Cambridge University Press. https://doi.org/10.1017/cbo9781139519328.011.

Hartigan, J. A., and M. A. Wong. 1979. Algorithm AS 136: A k-means clustering algorithm. Applied Statistics 28. JSTOR: 100. https://doi.org/10.2307/2346830.

Križ, Manuel, and Emmanuel Chemla. 2015. Two methods to find truth-value gaps and their application to the projection problem of homogeneity. Natural Language Semantics 23. Springer: 205–248. https://doi.org/10.1007/s11050-015-9114-z.

Lenth, Russell V. 2024. Emmeans: Estimated marginal means, aka least-squares means.

Strawson, P. F. 1964. Identifying reference and truth-values. Theoria 30. Wiley: 96–118. https://doi.org/10.1111/j.1755-2567.1964.tb00404.x.

Sun, Siqi, Karen M. Schmidt, and Teague R. Henry. 2025. Don’t let your Likert scales grow up to be visual analog scales: Understanding the relationship between number of response categories and measurement error. Ms; arXiv. https://doi.org/10.48550/ARXIV.2502.02846.

Thalmann, Maik. in progress. Not so squeamish, Your majesty: The presuppositions of singular definites between existence and uniqueness.

Thalmann, Maik, and Andrea Matticchio. submitted. Kleene’ing up your attitude: Against universal presupposition projection.

Thalmann, Maik, and Andrea Matticchio. 2024. On being certain that presuppositions don’t project universally. In Proceedings of the Amsterdam Colloquium, 378–385.

Appendix

German sample item

Alle Frauen sind (nicht) nochmal Kanu gefahren.
1. true:
  Maria ist in der Vergangenheit Kanu gefahren und vor Kurzem.
  Johanna ist in der Vergangenheit Kanu gefahren und vor Kurzem.
  Susanne ist in der Vergangenheit Kanu gefahren und vor Kurzem.
2. false:
  Maria ist in der Vergangenheit Kanu gefahren und vor Kurzem.
  Johanna ist in der Vergangenheit Kanu gefahren und vor Kurzem.
  Susanne ist in der Vergangenheit Kanu gefahren aber vor Kurzem nicht.
3. undefined:
  Maria ist in der Vergangenheit nicht Kanu gefahren aber vor Kurzem.
  Johanna ist in der Vergangenheit nicht Kanu gefahren und vor Kurzem nicht.
  Susanne ist in der Vergangenheit Kanu gefahren und vor Kurzem.
4. critical:
  Maria ist in der Vergangenheit nicht Kanu gefahren aber vor Kurzem.
  Johanna ist in der Vergangenheit Kanu gefahren und vor Kurzem.
  Susanne ist in der Vergangenheit Kanu gefahren aber vor Kurzem nicht.

Determining the optimal number of clusters

The question is how to determine the optimal number of clusters (K) for K-means clustering. Two common methods are the Elbow Method and the Silhouette Method. Both of them will be shown graphically.

The Elbow Method works by running the clustering algorithm repeatedly, each time with a different number of clusters ranging from \(1\) to \(16\), and recording how tightly packed the resulting clusters are. This tightness is measured by the total within-cluster sum of squares—essentially the sum of all the distances between each point and the centre of its assigned cluster. When plotted, this value drops steeply as the number of clusters increases from one, but at some point the curve flattens out. The optimal number of clusters is identified at the “elbow” of this curve: the point where adding more clusters stops producing meaningfully tighter groupings and therefore offers diminishing returns.

The Silhouette Method takes a different approach by asking, for each data point, how well it fits into its assigned cluster compared to how well it would fit into the next best alternative cluster. This produces a score between \(-1\) and \(1\) for every point, where higher values indicate a point is well-matched to its own cluster and poorly matched to neighbouring ones. The average of these scores across all points is then computed for each candidate number of clusters, and the number of clusters that produces the highest average silhouette width is considered optimal. Unlike the Elbow Method, which focuses purely on internal compactness, the Silhouette Method also rewards solutions where clusters are clearly separated from one another, making it a useful complementary check.

Code

mds_d <- cmdscale(dist(d_pivoted[, 2:length(d_pivoted)]))

# elbow Method for k-means clustering, compute WCSS for K = 1 to 16
wcss <- numeric()
for (k in 1:16) {
    set.seed(17)
    km <- kmeans(mds_d, centers = k, nstart = 25)
    wcss[k] <- km$tot.withinss
}
elbow_df <- data.frame(
    K = 1:16,
    WCSS = wcss
)

# silhouette method, compute average silhouette width for K = 2 to 16
avg_sil <- numeric()
for (k in 2:16) {
    set.seed(17)
    km <- kmeans(mds_d, centers = k, nstart = 25)
    sil <- silhouette(km$cluster, dist(mds_d))
    avg_sil[k] <- mean(sil[, 3])
}
sil_df <- data.frame(
    K = 2:16,
    AvgSilhouette = avg_sil[2:16]
) %>%
    add_row(K = 1, AvgSilhouette = NA)


p_elbow <- ggplot(elbow_df, aes(x = K, y = WCSS)) +
    geom_line(color = colors[1], linewidth = 1) +
    geom_point(color = colors[1], size = 3) +
    geom_vline(xintercept = 3, color = colors[2], lty = "dashed") +
    labs(
        title = "*Elbow Method* for determining optimal *K*",
        subtitle = "The 'elbow' where WCSS stops decreasing sharply identifies the number of clusters.",
        x = "Number of clusters (*K*)",
        y = "Within-cluster sums of squares"
    ) +
    scale_x_continuous(breaks = 1:16)

p_silhouette <- ggplot(sil_df, aes(x = K, y = AvgSilhouette)) +
    geom_line(color = colors[1], linewidth = 1) +
    geom_point(color = colors[1], size = 3) +
    geom_vline(xintercept = 2, color = colors[2], lty = "dashed") +
    labs(
        title = "*Silhouette Method* for determining optimal *K*",
        subtitle = "Average silhouette width near 1 = good clustering; near 0 = overlapping clusters.<br>The peak in the silhouette indicates the optimal number of clusters",
        x = "Number of clusters (*K*)",
        y = "Avg. silhouette width"
    ) +
    scale_y_continuous(limits = c(0.25, .75)) +
    scale_x_continuous(breaks = 1:16)


p_elbow / p_silhouette + plot_layout(axes = "collect")

Clustered by-participant plot

Here I show the descriptive results for an alternative clustering procedure that produces 2 instead of 3 clusters. As the plot indicates, the results are substantially similar, and as far as I can see, the decision does not meaningfully affect the conclusions I presented above.

Code

d_ac %>%
    arrange(cluster) %>%
    mutate(
        scenario = gsub("efined|ical", "", scenario),
        scenario = fct_relevel(scenario, "true", "false", "und"),
        negation = fct_relevel(negation, "without"),
        trigger = paste0("*", trigger, "*"),
        color = colors5[as.numeric(cluster)],
        cluster = fct_inorder(cluster),
        id = glue("<i style='color:{color}'>ID: {id} ({cluster})</i>"),
        id = fct_inorder(id)
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        pch = negation,
        color = negation,
        fill = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(
        fun.data = mean_cl_normal,
        geom = "errorbar",
        width = .1,
        alpha = .5
    ) +
    stat_summary(fun = mean, geom = "point", size = 2) +
    facet_wrap(~id, ncol = 6) +
    scale_y_continuous(breaks = c(-2, 0, 2)) +
    coord_cartesian(clip = "off", ylim = c(-2, 2)) +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        fill = "Negation",
        pch = "Negation"
    ) +
    guides(shape = guide_legend(override.aes = list(size = 3))) +
    theme_mt(base_size = 8, strip_text_size = 8) +
    theme(
        legend.margin = margin(-15, 0, 0, 0, "pt"),
        panel.spacing = unit(.25, "lines"),
        axis.text.x = element_markdown(family = "Cascadia Code", size = 6)
    )

Results with two clusters

Code

n_cluster <- 2
mds_d <- cmdscale(dist(d_pivoted[, 2:length(d_pivoted)]))
# Form clusters using K-means clustering
# changing the seed only affects the labelling of the clusters, not the sorting
set.seed(17)
kmeans_clusters_2 <- kmeans(
    mds_d,
    iter.max = 5000,
    algorithm = "Hartigan-Wong",
    nstart = length(unique(d$id)),
    centers = n_cluster
)$cluster

mds_df_2 <- as.data.frame(mds_d) %>%
    mutate(
        cluster = as.factor(kmeans_clusters_2),
        id = d_pivoted$id
    )
d_ac_2 <- d %>%
    left_join(
        mds_df_2 %>%
            select(id, cluster),
        by = "id"
    )

p_cluster_2 <- mds_df_2 %>%
    ggplot(aes(
        x = V1,
        y = V2,
        color = cluster,
        fill = cluster,
        shape = cluster
    )) +
    geom_point(size = 1.5) +
    geom_encircle(
        alpha = .2,
        s_shape = 1,
        expand = 0,
        show.legend = FALSE,
        color = NA
    ) +
    geom_text_repel(
        data = filter(mds_df_2, cluster == "1"),
        aes(label = id),
        box.padding = 0.5,
        direction = "x",
        seed = 1234,
        size = 3,
        min.segment.length = 1,
        max.iter = 50000,
        max.time = 2,
        show.legend = FALSE
    ) +
    scale_shape_manual(values = c(16, 15, 17)) +
    labs(
        x = "Multidimensional scaling (MDS) for participant clusters; dimension 1",
        y = "MDS dimension 2",
        color = "Cluster",
        fill = "Cluster",
        shape = "Cluster"
    ) +
    scale_color_manual(values = colors5) +
    scale_fill_manual(values = colors5) +
    guides(shape = guide_legend(override.aes = list(size = 3))) +
    theme(legend.margin = margin(-10, 0, 0, 0, "pt"))


p_judg_2 <- d_ac_2 %>%
    arrange(cluster) %>%
    mutate(
        scenario = fct_relevel(scenario, "true", "false", "undef"),
        negation = fct_relevel(negation, "without"),
        trigger = paste0("*", trigger, "*"),
        color = colors5[as.numeric(cluster)],
        cluster = glue("<i style='color:{color}'>**Cluster {cluster}**</i>"),
        cluster = fct_inorder(cluster)
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        pch = negation,
        color = negation,
        fill = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = TRUE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(fun.data = mean_cl_normal, geom = "errorbar", width = .1) +
    stat_summary(fun = mean, geom = "point", size = 2) +
    facet_wrap(trigger ~ cluster, ncol = n_cluster) +
    coord_cartesian(clip = "off") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        fill = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    guides(shape = guide_legend(override.aes = list(size = 3))) +
    theme(
        legend.margin = margin(-25, 0, 0, 0, "pt"),
        axis.text.x = element_markdown(family = "Cascadia Code", size = 7)
    )

p_cluster_2

Code

p_judg_2

Projection from attitude predicates

This experiment is reported in Thalmann and Matticchio (2024);Thalmann and Matticchio (submitted).

The experiment was set up as follows:

\(2 \times 4\) design
1. negation: without vs. with matrix negation (within items and participants)
2. scenario: true vs. false vs. undefined vs. critical (within items and participants)
3. Presupposition trigger as a pseudo-factor: again vs. stop (between items, 24 each)

Peter is (not) certain that Jan canoed again.
1. true Peter: “I’m certain that Jan canoed last time, and I’m certain that Jan canoed this time.”
2. false Peter: “I’m certain that Jan canoed last time, but I have no idea if this time, Jan canoed or not.”
3. undefined Peter: “I have no idea if last time, Jan canoed or not, but I am certain that Jan canoed this time.”
4. critical Peter: “I have no idea if last time, Jan canoed or not, and I have no idea if this time, Jan canoed or not.”

Code

read_csv(
    here("data", "believe-projection.csv"),
    show_col_types = FALSE
) %>%
    filter(sub_exp == 1) %>%
    mutate(
        scenario = gsub("undef", "undef.", scenario),
        scenario = fct_relevel(scenario, "true", "false", "undef."),
        negation = forcats::fct_recode(negation, without = "pos", with = "neg"),
        negation = fct_relevel(negation, "without", "with"),
        trigger = paste0("*", trigger, "*"),
        trigger = fct_relevel(trigger, "*again*"),
        judgment = judgment - 50,
        judgment = judgment / 25
    ) %>%
    ggplot(aes(
        x = scenario,
        y = judgment,
        color = negation,
        pch = negation,
        group = negation
    )) +
    geom_half_violin(
        data = . %>% filter(negation == "without"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[1],
        bw = .4,
        color = "transparent",
        trim = FALSE,
        scale = "count"
    ) +
    geom_half_violin(
        data = . %>% filter(negation == "with"),
        mapping = aes(group = scenario),
        adjust = .7,
        alpha = .2,
        fill = colors[2],
        bw = .4,
        color = "transparent",
        trim = FALSE,
        scale = "count",
        side = "R"
    ) +
    stat_summary(fun = mean, geom = "line") +
    stat_summary(
        fun.data = mean_cl_normal,
        width = .15,
        alpha = .7,
        geom = "errorbar"
    ) +
    stat_summary(fun = mean, geom = "point", size = 4) +
    guides(group = "none") +
    labs(
        x = "Scenario",
        y = "Judgment \u00B1CI<sub>95</sub>",
        color = "Negation",
        pch = "Negation"
    ) +
    scale_y_continuous(limits = c(-2, 2)) +
    facet_wrap(~trigger) +
    theme(axis.text.x = element_markdown(family = "Cascadia Code"))

Footnotes

This strategy is unlikely to be local accommodation since both the items with soft trigger stop and the items with the hard trigger again are affected. Ordinarily, we expect local accommation to be available only for the presuppositions of soft triggers.↩︎