How CBAS Works

An interactive walkthrough of the Choice-Wide Behavioral Association Study algorithm.

1 Counting Sequences

Imagine you're watching an animal make choices in a maze. At each decision point it goes Left (L) or Right (R). Over hundreds of trials you get a long stream of choices.

CBAS asks: do certain patterns appear more often in one group of subjects than another? To answer that, it first needs to count every possible short pattern in each subject's stream.

3

The sliding window moves across the stream and extracts every subsequence of the chosen length. Each unique pattern gets a count. With 2 symbols and length 3, there are 8 possible patterns (LLL, LLR, LRL, ...). With longer sequences or more symbols, the space grows exponentially.

This happens for every subject independently. At the end, you have a matrix where each row is a subject and each column is a sequence, and the values are how many times that subject produced that pattern.

The criterion parameter limits how much of the stream to use. Only windows starting at positions 0 through criterion (inclusive) are counted. This ensures all subjects contribute the same amount of data regardless of how long their sessions were.

2 Measuring Group Differences

Now you have a count for every sequence for every subject. You split subjects into two groups (e.g. control vs lesion) and ask: for each sequence, do the groups use it at different rates?

The intuition

Suppose 20 control rats produce the pattern "LRL" an average of 12 times each, while 20 lesion rats produce it only 8 times each. The groups differ by 4. But is that a meaningful difference, or just noise?

It depends on how much variation there is within each group. If individual rats range from 2 to 25 uses, a difference of 4 isn't impressive. If they range from 10 to 14, a difference of 4 is huge. The test statistic captures this by dividing the group difference by a measure of the spread.

Building the t-statistic step by step

1. Compute the mean for each group

Add up the counts for all subjects in a group and divide by how many subjects there are. This gives you the average usage rate for this sequence in each group.

2. Compute the standard error of the mean (SEM) for each group

The SEM tells you how precisely you know the group mean. It's the standard deviation divided by the square root of the sample size. A larger group or less variable group gives you a smaller SEM (more precision).

3. Combine the SEMs

Because you're comparing two uncertain quantities (both group means have noise), you combine their uncertainties: √(SEM₀² + SEM₁²). This is like error propagation in physics.

4. Divide the difference by the combined uncertainty

The final statistic: how many "combined standard errors" apart are the groups?

t = (meangroup 0 − meangroup 1) / √(SEM0² + SEM1²) A t of 1 means the groups are 1 combined-SEM apart. A t of 4 means they're 4 apart. Bigger = more evidence of a real difference.
t = 3.2

Try dragging both sliders. When you increase separation, t goes up. When you increase spread (noise within each group), t goes down. The statistic captures the signal-to-noise ratio.

Why two one-sided tests?

For each sequence, CBAS only stores the t-statistic in the direction that was actually observed. If group 0 used it more, t goes in the "positive" slot. If group 1 used it more, the absolute value of t goes in the "negative" slot. The other slot is left empty (NaN).

This doubles the number of hypotheses (from S sequences to 2S tests) but lets the algorithm handle directionality cleanly. A sequence is significant if it differs in either direction, and the output tells you which.

Think of it like measuring wind. You record whether wind blows east or west, and how fast. You don't mix eastward and westward measurements into a single number because they'd cancel out. Instead you track each direction separately.

3 Building the Null Distribution

You have a large t-value for some sequence. But is it large enough? When you test thousands of sequences simultaneously, some will have big t-values purely by chance. You need to know what "big by chance" looks like.

The shuffling trick

The bootstrap answers: "If there were truly no difference between the groups, what would the test statistics look like?" It finds out by erasing the group labels, pooling everyone together, randomly assigning subjects to fake groups, and recomputing the statistics.

Each time you shuffle and recompute, you get one "null" draw for every sequence simultaneously. Do this 10,000 times and you get a matrix of null values.

The bootstrap matrix

This matrix is the central data structure for everything that follows. It has 10,000 rows (one per resample) and as many columns as there are sequences being tested. Each cell contains "what this sequence's t-value would have been under random group assignment number B."

The columns are sequences (the same ones you computed real t-values for). The rows are resamples (each one a different random shuffling of group labels). When we talk about "a row" later, we mean one complete set of null t-values from one shuffling. When we talk about "a column," we mean one sequence's null distribution across all shuffles.

From the matrix to a p-value for one sequence

For a single sequence, you could get its p-value by looking down its column and asking "what fraction of null values are bigger than my observed value?" But that ignores the multiple testing problem. The step-down procedure (step 4) uses the matrix differently, looking across rows rather than down columns.

0 / 10,000

This histogram shows the distribution of row maximums (the biggest null t-value in each row). The red dashed line is your real observed statistic. The fraction of row-maximums that exceed it is the adjusted p-value. This is much stricter than looking at just one column because the row-maximum represents "what the most extreme random noise across all sequences would look like."

All sequences share the same shuffled groups on each draw. This preserves correlations between sequences (e.g. if LLL goes up, LLLL probably goes up too). Preserving this structure matters for the step-down correction next.

4 Correcting for Multiple Comparisons

Here's the core problem. If you test 2,000 sequences and use a p-value threshold of 0.05, you'd expect about 100 false positives even if nothing is real. You need a way to adjust for the fact that you ran so many tests.

The naive approach (Bonferroni)

The simplest fix: divide your significance threshold by the number of tests. With 2,000 tests and alpha = 0.05, you'd need p < 0.000025 to declare significance. This works but is extremely conservative. Real effects with moderate strength get missed because the bar is so high.

The step-down idea

Romano-Wolf step-down is a smarter approach. It works in rounds, peeling off the most significant sequences one at a time, and each round uses a slightly easier threshold.

Here's the logic:

  1. Sort all your observed test statistics from largest to smallest.
  2. Start with the biggest one. Ask: "Is this larger than what I'd expect from pure chance?" To answer, look at each row in the bootstrap matrix. In each row, find the maximum value across all columns (sequences). You now have 10,000 "row maximums." If your observed statistic beats most of them, it's significant.
  3. Remove it and move on. Now look at the second-largest statistic. But this time, the column for the sequence you just declared significant is ignored when computing row maximums. So the maximum in each row is computed over the remaining columns only. This makes the comparison slightly easier.
  4. Repeat until a statistic fails to beat its threshold. Everything from that point down also fails (monotonicity).

Where does the threshold come from?

The threshold at each step is calculated from the bootstrap matrix:

  1. Take all 10,000 rows. For each row, find the maximum across the remaining columns (sequences not yet rejected).
  2. Sort those 10,000 row-maximums from smallest to largest.
  3. Pick the value at position 10,000 × (1 − α). With α = 0.5, that's the median (position 5,000).

Your observed statistic must beat this quantile to be declared significant at this step.

Why does removing a column lower the threshold? Because the maximum over 1,999 columns tends to be smaller than the maximum over 2,000 columns. Every time you peel off a rejected sequence's column, each row's maximum is computed over a smaller set, pulling the whole distribution of row-maximums downward. The quantile of that distribution drops too.
Imagine a spelling bee. In each round, the best speller is tested against the hardest word. If they get it right, they're declared "exceptional" and removed. The next round uses slightly easier words because the exceptional speller is gone. Each round, the threshold adapts to who's left.

The blue bars are observed test statistics (sorted largest to smallest). The red dashed line is the threshold at each step, computed as the median row-maximum over the remaining columns. Watch how it drops as rejected sequences (pink) have their columns removed. Grey bars failed to beat the threshold.

Why this is more powerful than Bonferroni

Bonferroni sets one fixed threshold for all tests. Step-down gives a different (lower) threshold to each successive sequence. The strongest signals clear the high bar, and once they're gone, moderate signals have an easier bar to clear. This matters a lot when you have hundreds of truly significant sequences mixed in with thousands of null ones.

Direction-conditional removal

One subtlety: when a sequence's column is "removed" from a bootstrap row, it only actually disappears if the bootstrap happened to go the same direction as the real observation. If the bootstrap went the opposite direction for that sequence in that row, the value stays. This is a conservative choice that prevents the procedure from being too generous.

5 Allowing Some False Positives (k-FWER)

The step-down procedure as described above controls the "family-wise error rate" at k=1. That means: the probability of even one false positive is kept below alpha. This is very strict.

Why strictness hurts

If you truly have 800 significant sequences, insisting on zero false positives forces an extremely high bar. You might only recover 400 of them. Is it worth missing 400 real effects just to avoid a single false one?

The k-FWER relaxation

The fix is simple: instead of computing the row maximum (the single biggest null value in each row), compute the k-th largest value. With k=1, you get the maximum. With k=5, you get the 5th largest. A higher k means a lower threshold, which means more sequences pass.

The tradeoff: you're now tolerating up to k false positives among your rejections. The threshold is lower because you're essentially saying "I'm okay if up to k of my rejections are wrong."

k = 1

The bars are the same observed test statistics from before. The threshold (red line) is now the quantile of the k-th largest per row instead of the row maximum. Drag k higher and watch the threshold drop and more bars turn pink (significant). The number of rejections increases as you relax the standard.

Finding the right k automatically

You don't pick k by hand. The algorithm finds the k where the expected fraction of false positives among your rejections is controlled at level γ (default 0.05, meaning at most 5% false).

The logic is iterative:

  1. Start with k=1. Run the full step-down. Count how many rejections R you get.
  2. Check the balance. If R rejections and at most k of them are false, the false discovery proportion is k/R. Is k/R ≤ γ? In other words, is R ≤ k/γ?
  3. If R is too large (meaning k/R < γ, so you have too many rejections relative to your tolerance), increase k. Specifically, set the new k to the smallest value where k/R ≥ γ would hold.
  4. Repeat. A larger k gives a lower threshold, which gives more rejections R, which might require an even larger k. The process converges when R and k come into balance.
Think of it like adjusting the sensitivity on a metal detector at the beach. At maximum sensitivity (k=1), you only dig when the signal is overwhelming, so you never dig up bottle caps but you also walk past coins buried a bit deeper. Turning down the sensitivity (higher k) means you dig more often and find more coins, but also start pulling up some junk. You keep adjusting until the ratio of junk to total finds stabilizes at your tolerance level (5%).
γ = 0.05

The left panel shows each iteration: the k tried, how many rejections R resulted, the ratio k/R, and whether k/R ≥ γ. The right panel plots the trajectory. Each dot is one iteration, moving rightward as k increases. The green shaded region is where k/R ≥ γ (the acceptable zone). The algorithm stops the moment a dot lands in the green zone, because at that k the false discovery proportion is under control. Try adjusting γ to see how a more tolerant threshold converges at a different k.

What γ means in practice

With γ = 0.05 and 1,605 rejections, you're accepting that about 81 of them (5%) might be false. This is the k that the fly dataset converges to. In exchange for tolerating those ~81 potential false positives, you recover hundreds more true positives that would have been missed with k=1.

6 The Full Picture

CBAS takes behavioral streams, tests every short pattern for group differences, and returns adjusted p-values (called ζ-values) that account for the massive multiple testing burden.

The power of the approach comes from three things working together:

A sequence is declared significant if its ζ-value is below α (default 0.5). The final output includes the direction of the effect, letting you say not just "these groups differ on this pattern" but "group 0 uses this pattern more than group 1."

← Technical algorithm description