Convolutional Regression Tsetlin Machine — reading a number off an image¶
Paper. K. Darshana Abeyrathna, Ole-Christoffer Granmo, Morten Goodwin, Convolutional Regression Tsetlin Machine: An Interpretable Approach to Convolutional Regression, ICMLT 2021. doi:10.1145/3468891.3468901 (PDF in
source_documents/papers/03_segmentation/3468891.3468901.pdf)
This is one of only two papers in the repository's 03_segmentation folder, and it is not
a segmentation paper. It earns its place there because it is the closest published mechanism
for making a Tsetlin machine emit something other than a class label: a continuous value,
computed from patches. If you want a TM to produce a dense output — a mask, a depth map, a
per-pixel score — this is the machinery you would build on.
The C-RTM is an ensemble of two earlier ideas:
| from | what it contributes |
|---|---|
| Convolutional TM (Granmo et al., 2019) | clauses act as filters that slide over the image; a clause is True if it matches at least one patch |
| Regression TM (Abeyrathna et al., 2020) | clause polarity is dropped, every clause is additive, and the clause count is rescaled into a real number |
torchtsetlin ships this model as ConvRegressionTsetlinMachine, so the whole paper is
three lines of setup. What this notebook does with it:
- rebuild the paper's artificial datasets exactly (Figure 3 and Table 1);
- check the library's output layer and feedback rule against Equations 5 and 6 algebraically;
- reproduce the paper's headline claim — zero MAE on noise-free data — and verify the Table 5 arithmetic that produces it. With the current library this needs fewer clauses and longer training than the paper's settings; section 4 says why;
- decode the clauses into statements about named pixels, and show the machine recovered the data-generating mask — while finding a different decomposition from the paper's;
- be honest about how narrow that exact-recovery regime is: measure it over restarts and mask placements, and give a selection rule that needs no test labels;
- run the noisy-data comparison against two CNNs (Tables 6 and 7) — where the paper's C-RTM result does not reproduce, which section 9 documents rather than hides.
Runtime: about 6 minutes on a GPU (sections 7 and 8 fit 15 models between them).
Written for a GPU (about 6 minutes on an RTX 3080; the longest part is the restart grid in section 7). In Colab choose Runtime → Change runtime type → T4 GPU; expect it to take noticeably longer there.
# Setup: on Google Colab this installs torchtsetlin from PyPI (torch, torchvision, scikit-learn,
# matplotlib and pandas are already there). Elsewhere: `pip install torchtsetlin`, or
# `pip install -e ".[dev]"` in a checkout of the repository.
import sys
IN_COLAB = "google.colab" in sys.modules
if IN_COLAB:
%pip install -q "torchtsetlin>=0.7"
import itertools
import time
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import torch
import torch.nn as nn
import torchtsetlin as tt
SEED = 0
tt.seed_everything(SEED)
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
print("torchtsetlin", tt.__version__, "| torch", torch.__version__, "| device", device)
1. The data — the paper's Figure 3, exactly¶
The paper builds its datasets from a weight-mask. An image is X × Y × 1 with uniformly
random binary pixels. A mask marks k important pixels and gives them the place values
2^(k-1), …, 2^1, 2^0. The regression target is
$$\hat{y} \;=\; 100 \times \sum_{i \in \text{important}} 2^{\,e_i}\, x_i .$$
So for k = 4 the target takes the 16 values 0, 100, 200, …, 1500, and — this is the point
of the construction — the input–output mapping is known exactly, which is what lets the
paper claim optimality rather than just good accuracy.
Table 1 of the paper enumerates the 72 datasets: image sizes 3×3, 4×4, 5×5 × k ∈ {2,3,4}
× 4 random mask placements × {noise-free, noisy}.
def make_weighted_mask_dataset(n, size=4, k=4, seed=0, noise=0.0, device=device):
"""Paper Fig. 3: random binary images + a k-pixel weight-mask with place values 2^(k-1)..2^0.
Returns ``x`` (n,1,size,size), the clean target ``y``, the (possibly noisy) training
target ``y_noisy``, and the mask itself.
"""
g = torch.Generator().manual_seed(seed)
coords = torch.randperm(size * size, generator=g)[:k]
place = torch.tensor([2.0 ** (k - 1 - i) for i in range(k)])
mask = torch.zeros(size * size)
mask[coords] = place
mask = mask.view(size, size)
x = (torch.rand(n, 1, size, size, generator=g) < 0.5).float()
y = 100.0 * (x[:, 0] * mask).sum(dim=(1, 2))
y_noisy = y + (torch.rand(n, generator=g) * 2 - 1) * noise if noise else y.clone()
return x.to(device), y.to(device), y_noisy.to(device), mask
SIZE, K = 4, 4
x, y, _, mask = make_weighted_mask_dataset(4000, SIZE, K, seed=1)
x_train, y_train, x_test, y_test = x[:3000], y[:3000], x[3000:], y[3000:]
Y_MAX = 100.0 * (2 ** K - 1)
print("mask (place values):")
print(mask.int().numpy())
print(f"\nimages {tuple(x.shape)} targets in {{{', '.join(str(int(v)) for v in sorted(set(y.tolist()))[:4])}, ..., {int(Y_MAX)}}}")
print("distinct targets:", len(set(y.tolist())), "= 2^k")
# The paper's Figure 3, drawn for one example: image x mask -> resulting image -> output.
i = int((y[:200] == 1300).nonzero()[0]) if (y[:200] == 1300).any() else 0
fig, axes = plt.subplots(1, 3, figsize=(9.5, 3.0))
for ax, img, title in zip(
axes,
[x[i, 0].cpu(), mask, x[i, 0].cpu() * mask],
["image $x$ (random bits)", "weight-mask", "resulting image"],
):
ax.imshow(img, cmap="Blues", vmin=0, vmax=8)
for (r, c), v in np.ndenumerate(img.numpy()):
ax.text(c, r, f"{v:g}", ha="center", va="center",
color="white" if v > 4 else "black", fontsize=11)
ax.set_title(title, fontsize=10)
ax.set_xticks([]); ax.set_yticks([])
fig.suptitle(f"sum of the resulting image x 100 = y = {y[i]:.0f}", y=0.02, fontsize=11)
plt.tight_layout()
plt.show()
2. From CTM to C-RTM — the three changes¶
2.1 Polarity is dropped¶
A classification TM splits clauses into positive (vote for) and negative (vote against)
halves, so the vote sum lives in [-T, T]. The C-RTM removes polarity: every clause is
additive, and the vote sum lives in [0, T].
2.2 The output layer (paper Eq. 5)¶
$$y_o \;=\; \frac{\sum_{j=1}^{m} c_j(\hat{X}_o)}{T} \times \hat{y}_{\max}$$
torchtsetlin's RegressionTsetlinMachine implements the slightly more general
y = y_min + (v/T)·(y_max − y_min). With y_range=(0, ŷmax) the two are identical, and T
is the output resolution: the model can only emit multiples of ŷmax / T.
2.3 The feedback rule (paper Eq. 6)¶
Type I if
y_o < ŷ_o· Type II ify_o > ŷ_o, with probabilityK·|y_o − ŷ_o| / ŷmax.
Type I forces clauses towards firing, Type II forces them away — so the output is pushed
up when it is too low and down when it is too high. The library selects each clause with
probability ((v* − v) / T)², where v* is the target expressed in vote units. Substituting
y = (v/T)·ŷmax shows what that is in the paper's terms:
$$\frac{|y_o - \hat{y}_o|}{\hat{y}_{\max}} \;=\; \frac{\bigl|\tfrac{v}{T}\hat{y}_{\max} - \tfrac{v^\ast}{T}\hat{y}_{\max}\bigr|}{\hat{y}_{\max}} \;=\; \frac{|v - v^\ast|}{T},$$
so the library's probability is the square of the paper's Eq. 6 with K = 1. The square is
the rule of the reference implementations (pyTsetlinMachine, pyTsetlinMachineParallel and tmu),
which torchtsetlin follows since version 0.5. The two agree on the direction of every
feedback and at the extremes (no error, or an error of the whole range), but not in between:
an output that is off by one 100-unit step gets feedback with probability 1/15 ≈ 6.7 % under
Eq. 6 and 1/225 ≈ 0.44 % under the squared rule. Small errors are corrected far more slowly,
which matters in section 4.
Let's verify the output layer numerically rather than take it on faith.
probe = tt.ConvRegressionTsetlinMachine(
n_clauses=15, T=15, s=3.0, patch_size=(2, 2), stride=1,
position_encoding=True, y_range=(0.0, Y_MAX), input_shape=(1, SIZE, SIZE),
).to(device)
probe(x_test[:1]) # initialise lazily
v = torch.tensor([[0.0], [1.0], [7.0], [13.0], [15.0]], device=device)
paper = (v.squeeze(1) / probe.T) * Y_MAX # Eq. 5
library = probe.votes_to_targets(v) # torchtsetlin
print(pd.DataFrame({"vote sum v": v.squeeze(1).cpu(), "paper Eq. 5": paper.cpu(),
"votes_to_targets": library.cpu()}).to_string(index=False))
assert torch.allclose(paper, library)
print(f"\noutput resolution = y_max / T = {Y_MAX} / {probe.T:.0f} = {Y_MAX / probe.T:g}"
f" <- exactly the 100-unit grid the data lives on")
3. Position encoding — the paper's Table 4¶
A convolutional clause sees only a W × W patch, so on its own it cannot tell where it
matched. The paper's answer (Section 3.1, Table 4) is to append thermometer-encoded filter
coordinates to the patch literals, "to make the clause location-aware". Without them a
clause that means "the top-left pixel of my patch is 1" would fire wherever any such pixel
exists, which on a 4×4 image is almost always.
torchtsetlin builds the same bits in _ConvMixin._position_bits: P_y − 1 bits [y > k]
and P_x − 1 bits [x > k]. The paper tabulates ≤ thresholds instead of >; the two are
complements of each other and carry identical information.
Py, Px = probe._grid(SIZE, SIZE)
bits = probe._position_bits(Py, Px, torch.device("cpu")).int()
rows = [{"location (x, y)": f"({b % Px}, {b // Px})",
**{f"y>{k}": int(bits[b, k]) for k in range(Py - 1)},
**{f"x>{k}": int(bits[b, Py - 1 + k]) for k in range(Px - 1)}}
for b in range(Py * Px)]
print(f"patch grid {Py} x {Px} = {Py * Px} filter positions "
f"(paper: B = BX x BY = {Py} x {Px})\n")
print(pd.DataFrame(rows).to_string(index=False))
print(f"\npatch literals: {K} pixels + {(Py - 1) + (Px - 1)} position bits "
f"= {probe.n_features} features -> {2 * probe.n_features} literals")
4. Training the C-RTM¶
The paper picks the clause budget from the Regression-TM analysis: 2^k − 1 clauses
(3 for k=2, 7 for k=3, 15 for k=4) with T = m. The arithmetic behind that number is
Table 5: to emit 1300 the machine must fire clauses worth 8 + 4 + 1 place values, so it
needs 8 clauses that recognise the weight-8 pixel, 4 for the weight-4 pixel, 2 for weight-2
and 1 for weight-1 — 8 + 4 + 2 + 1 = 15 in total. The paper calls this weighting
by duplication.
torchtsetlin also offers the Integer-Weighted TM (weighted=True), where a clause carries
an integer weight instead of being copied. Both are worth running, because they behave very
differently here.
The duplication model below uses the paper's settings (15 clauses, T = 15). The
integer-weighted model keeps T = 15 — that is the output resolution, the 100-unit grid — but
uses 8 clauses, 400 epochs and batches of 16. With the paper's 15 weighted clauses (and the
100–150 epochs of batch 8 used here before) the current library did not reach the exact
solution in any of 12 restarts; the test MAE stayed around 75–230. What we tried and what
helped:
integer-weighted setting (T = 15, s = 3, batch 8 unless stated) |
restarts with test MAE 0 (4 masks × 3 seeds) |
|---|---|
| 15 clauses, 100 epochs, batch 8 | 0 / 12 |
15 clauses: 300 epochs; or, at 100 epochs, s = 2 or 5, no clause budget, or batch 1 |
0 / 4 each (one seed per mask) |
| 10 clauses, 100 epochs / 300 epochs | 0 / 12 · 6 / 12 |
| 8 clauses, 100 / 200 / 400 epochs | 1 / 12 · 5 / 12 · 8 / 12 |
| 8 clauses, 400 epochs, batch 16 (used below) | 8 / 12 |
Why fewer clauses and more epochs? Two rules changed in torchtsetlin 0.5 to match the
reference implementations: a clause weight can no longer drop below 1 (it used to be able to
reach 0), and the feedback probability became the squared error of section 2.3. A clause that
is not needed can therefore no longer be silenced by its weight — it has to learn never to fire
— and the squared rule makes the last one-step corrections rare. Our reading (not separately
tested) is that fewer spare clauses and more epochs simply give both effects the room they
need. An earlier version of this notebook, on torchtsetlin 0.1.x, found
the exact solution with the paper's 15 clauses on 9 of 12 restarts. In side experiments (not
part of this notebook) we checked that these two rules are the cause: putting both old rules
back into the current library restored exact recovery on 6 of 12 restarts with the old settings,
either one alone on none, and the 0.1.2 code with both new rules found it on none.
def fit_crtm(x_tr, y_tr, *, n_clauses, T, s=3.0, weighted=False, epochs=150, batch_size=8,
seed=0, patch_size=(2, 2), size=SIZE, y_max=Y_MAX, max_lit=5, pos=True):
tt.seed_everything(seed)
model = tt.ConvRegressionTsetlinMachine(
n_clauses=n_clauses, T=T, s=s, patch_size=patch_size, stride=1,
position_encoding=pos, y_range=(0.0, y_max), input_shape=(1, size, size),
max_included_literals=max_lit, weighted=weighted,
).to(x_tr.device)
for _ in range(epochs):
perm = torch.randperm(x_tr.shape[0], device=x_tr.device)
for i in range(0, len(perm), batch_size):
model.update(x_tr[perm[i:i + batch_size]], y_tr[perm[i:i + batch_size]])
model.eval()
return model
def mae(model, xs, ys):
return (model.predict(xs) - ys).abs().mean().item()
# Integer-weighted settings (see the text below): T = 15 keeps the 100-unit output grid, but
# only 8 clauses, trained longer and with batches of 16.
N_CLAUSES_W, EPOCHS_W, BATCH_W = 8, 400, 16
t0 = time.time()
by_duplication = fit_crtm(x_train, y_train, n_clauses=15, T=15, weighted=False, seed=0)
by_weight = fit_crtm(x_train, y_train, n_clauses=N_CLAUSES_W, T=15, weighted=True,
epochs=EPOCHS_W, batch_size=BATCH_W, seed=0)
print(f"({time.time() - t0:.0f}s)\n")
print(f"{'weighting':<26s} {'train MAE':>10s} {'test MAE':>10s}")
for name, m in (("by duplication (paper)", by_duplication), (f"integer weights ({N_CLAUSES_W} cl.)", by_weight)):
print(f"{name:<26s} {mae(m, x_train, y_train):>10.2f} {mae(m, x_test, y_test):>10.2f}")
On this dataset the integer-weighted model typically reaches an MAE of zero — the paper's headline claim, "The C-RTM obtains zero MAE for both training and testing on noise-free datasets" — though, as section 4's table shows, only with fewer clauses and longer training than the paper's, and section 7 shows it is not guaranteed for every restart. The duplication form, run with the paper's own budget, does not converge (its MAE stays in the hundreds).
That difference is worth understanding rather than papering over. Duplication needs fifteen
independent clauses to settle into the exact multiset {8 × A, 4 × B, 2 × C, 1 × D}, and every
clause receives the same error signal on every example — nothing in the feedback rule
differentiates them, so the split has to emerge from initialisation noise alone. An integer
weight reaches the same arithmetic with as few as four clauses and an explicit counter.
5. Reading the machine¶
This is the paper's real selling point over a CNN: the fifteen filters are printable. Table 2
of the paper shows patterns like (1 ⋆ ⋆ ⋆) — "this pixel must be 1, the rest are don't-care"
— and Table 3 shows them pinned to locations. Both fall straight out of clause_patch
(1 = must be on, -1 = must be off, 0 = don't care) and clause_region (the decoded
position literals).
from torchtsetlin import functional as Fnl
def clause_conditions(model, xs):
"""Decode every clause into the absolute-pixel condition it tests.
``clause_patch`` gives the within-patch pattern (``1`` = must be on, ``-1`` = must be off,
``0`` = don't care) and ``clause_region`` decodes the position literals into the patch
coordinates the clause is pinned to. Together they name real pixels of the image.
"""
lit = model._encode(model._prepare(xs))
B, P, L = lit.shape
fires = Fnl.clause_outputs(lit.reshape(B * P, L), model.include,
model.include_count, False).view(B, P, -1).any(dim=1).float()
kh, kw = model.patch_size
rows = []
for j in range(model.n_clauses_total):
w, rate = int(model.weights[j]), round(fires[:, j].mean().item(), 3)
if int(model.include_count[j]) == 0:
rows.append({"clause": j, "weight": w, "fires": rate, "condition": "(empty clause)"})
continue
reg = model.clause_region(j)
(ylo, yhi), (xlo, xhi) = reg["y"], reg["x"]
if ylo > yhi or xlo > xhi:
cond = "(contradictory — never fires)"
elif ylo != yhi or xlo != xhi:
cond = f"not pinned (y={ylo}..{yhi}, x={xlo}..{xhi})"
else:
pat = model.clause_patch(j).squeeze(0)
parts = [("" if int(pat[a, b]) > 0 else "NOT ") + f"({ylo + a},{xlo + b})"
for a in range(kh) for b in range(kw) if int(pat[a, b]) != 0]
cond = " AND ".join(parts) if parts else "(always true)"
rows.append({"clause": j, "weight": w, "fires": rate, "condition": cond})
return pd.DataFrame(rows)
err = (by_weight.predict(x_test) - y_test).abs()
print(f"test MAE {mae(by_weight, x_test, y_test):.4f} max absolute error {err.max():.2e} "
f"fraction within half a step (50): {err.lt(50).float().mean():.3f}\n")
print(clause_conditions(by_weight, x_test).to_string(index=False))
print("\nimportant pixels and their place values:",
{tuple(p): int(mask[tuple(p)]) for p in (mask > 0).nonzero().tolist()})
Every clause that fires is a readable statement about named pixels of the image — the paper's Table 2 and Table 3 combined, decoded automatically.
But the decomposition is usually not the idealised one the paper describes. The paper's
Table 5 assumes one pattern per important pixel (8 × (1 ⋆ ⋆ ⋆) + 4 × (⋆ 1 ⋆ ⋆) + …). Some
pixels do get a single clause (in our run the weight-1 and weight-2 pixels: (1,0) with
weight 1, (1,2) with weight 2). Others are covered by a case split: clauses whose
conditions are complementary in a second pixel, so exactly one of them fires whenever the
important pixel is on. In our run, for example:
pixel (1,1), place value 8: NOT (1,0) AND (1,1) -> 8
(1,0) AND (1,1) -> 6 + 2 (two clauses, same condition)
pixel (3,3), place value 4: NOT (3,2) AND (3,3) -> 4
NOT (2,2) AND (3,2) AND (3,3) -> 4
(2,2) AND (3,2) AND (3,3) -> 4
Each group contributes its place value precisely when the important pixel is on — arithmetically identical to a single clause of that weight, reached by a different route. The exact clauses and weights vary from run to run (and between GPUs), so expect a different decomposition of the same kind rather than this one. Note that every clause carries a weight of at least 1, so a clause that is not part of the solution cannot be switched off through its weight; it has to become a condition that never fires.
The lesson is that the exact solution is not unique, and the interpretability claim survives it: the clauses are still readable, they are just not in canonical form.
What the machine actually recovered¶
Clause weights alone therefore do not spell out the mask. The question the mask answers —
how much does this pixel contribute? — is answered by asking the model directly. Because
y = 100 Σ w_p x_p is additive, flipping one pixel from 0 to 1 must change the prediction by
exactly 100 w_p, whatever the other pixels are doing.
def sensitivity_map(model, xs):
"""Average change in prediction when each pixel is forced from 0 to 1."""
out = torch.zeros(SIZE, SIZE)
for r in range(SIZE):
for c in range(SIZE):
off, on = xs.clone(), xs.clone()
off[:, 0, r, c], on[:, 0, r, c] = 0.0, 1.0
out[r, c] = (model.predict(on) - model.predict(off)).mean().item()
return out
recovered = sensitivity_map(by_weight, x_test)
truth = 100.0 * mask
fig, axes = plt.subplots(1, 2, figsize=(7.2, 3.1))
for ax, img, title in zip(axes, [truth, recovered],
["true contribution (100 x mask)", "recovered by probing the model"]):
ax.imshow(img, cmap="Blues", vmin=0, vmax=800)
for (r, c), v in np.ndenumerate(img.numpy()):
ax.text(c, r, f"{v:.0f}", ha="center", va="center",
color="white" if v > 400 else "black", fontsize=12)
ax.set_title(title, fontsize=10); ax.set_xticks([]); ax.set_yticks([])
plt.tight_layout(); plt.show()
print("identical:", torch.equal(recovered.round(), truth))
When training reached the exact solution (test MAE 0 above), the machine did not merely fit
the function; it recovered the data-generating mask, place value for place value, including
the twelve pixels that contribute nothing. (A run that missed the exact solution shows up here
as identical: False.)
6. Reproducing Table 5¶
Table 5 of the paper lists, for every attainable output, which clauses have to fire. With
duplication that table has 16 rows of clause counts; with integer weights it collapses to
the binary expansion of y / 100. Either way the check is the same: enumerate all 2^k
input configurations of the important pixels and confirm the machine emits the right number.
important = (mask > 0).nonzero().tolist()
rows = []
for combo in itertools.product([0.0, 1.0], repeat=K):
img = torch.zeros(1, 1, SIZE, SIZE, device=device)
# fill the non-important pixels randomly to prove they really are don't-care
img[0, 0] = (torch.rand(SIZE, SIZE, device=device) < 0.5).float()
for (r, c), bit in zip(important, combo):
img[0, 0, r, c] = bit
target = 100.0 * sum(float(mask[r, c]) * b for (r, c), b in zip(important, combo))
pred = float(by_weight.predict(img))
lit = by_weight._encode(by_weight._prepare(img))
B, P, L = lit.shape
fired = Fnl.clause_outputs(lit.reshape(B * P, L), by_weight.include,
by_weight.include_count, False).view(B, P, -1).any(1)[0]
w = by_weight.weights[fired].tolist()
rows.append({"important pixels": "".join(str(int(b)) for b in combo),
"target": int(target), "predicted": round(pred, 1),
"clause weights firing": "+".join(map(str, sorted(w, reverse=True))) or "none",
"sum": int(sum(w))})
tab5 = pd.DataFrame(rows)
print(tab5.to_string(index=False))
print(f"\nall {len(tab5)} outputs exact:", bool((tab5['target'] == tab5['predicted']).all()))
7. How narrow is the exact-recovery regime?¶
The paper flags this itself, in a passage that is easy to skim past:
"We observe that the C-RTM behaves similar to RTM when T (equal to the number of clauses) is not a multiplier of the optimum required clauses. […] if there are five clauses for the pattern (⋆ 1 ⋆ ⋆) after training, [the] input image […] activates 14 clauses, producing an incorrect output that is 1400."
In other words the machine has to land on exactly the right clause configuration; a neighbouring one is not approximately right, it is off by a whole 100-unit step. Let's measure how often that happens across mask placements and initialisations, instead of reporting one good run.
t0 = time.time()
grid = []
for mask_seed in range(4):
xs, ys, _, _ = make_weighted_mask_dataset(4000, SIZE, K, seed=mask_seed)
xtr, ytr, xte, yte = xs[:3000], ys[:3000], xs[3000:], ys[3000:]
for tm_seed in range(3):
m = fit_crtm(xtr, ytr, n_clauses=N_CLAUSES_W, T=15, weighted=True,
epochs=EPOCHS_W, batch_size=BATCH_W, seed=tm_seed)
grid.append({"mask": mask_seed, "init": tm_seed,
"train MAE": round(mae(m, xtr, ytr), 3),
"test MAE": round(mae(m, xte, yte), 3)})
grid = pd.DataFrame(grid)
print(f"({time.time() - t0:.0f}s)\n")
print(grid.pivot(index="mask", columns="init", values="test MAE").to_string())
exact = (grid["test MAE"] < 0.01).sum()
print(f"\nruns reaching the exact decomposition (test MAE < 0.01): {exact}/{len(grid)}")
So the paper's "zero MAE" is real but it is a basin you fall into or miss, not a limit you converge to. In our run 8 of the 12 restarts found it — every restart on masks 0 and 1 and two of three on mask 2, but none on mask 3, whose four important pixels all sit inside one 2×2 window; the misses are off by roughly 50–170 MAE, i.e. whole 100-unit steps on part of the examples. Mask 3 is not out of reach in principle (the 0.1.x library found it), but none of the current-library settings in section 4's table did, including six restarts at 8 clauses and 400 epochs. On a different GPU the individual numbers will differ; the pattern — most restarts exact, a mask that resists — is what to expect.
Fortunately you do not need test labels to tell the outcomes apart: a run that found the exact decomposition has a training MAE of 0 as well, and a run that missed it does not. Restarting a handful of times and keeping the model with the lowest training error is therefore a legitimate recipe — it finds the exact solution whenever one of the restarts did, and reports honestly (through the training MAE) when none did.
print(grid.assign(**{"found exact": grid["test MAE"] < 0.01})
.groupby("found exact")[["train MAE", "test MAE"]].agg(["min", "max"]).to_string())
selected = grid.loc[grid.groupby("mask")["train MAE"].idxmin()]
print("\nper mask, selecting the restart with the lowest TRAINING MAE:")
print(selected.to_string(index=False))
print(f"\n-> exact on {(selected['test MAE'] < 0.01).sum()}/{len(selected)} masks")
8. Ablation — what the position bits are worth¶
Section 3.1 of the paper argues the location literals are what stop a clause from firing at
the wrong place ("why the clause which recognizes pattern 1 does not output 1 at picture
location (2, 0) […] these clauses also learn the correct location of the pattern").
torchtsetlin makes position_encoding a constructor flag, so the claim is one run away.
# Mask 1 of the section 7 grid is this notebook's main dataset (seed=1), trained with
# position_encoding=True and the same seeds, so only the runs without position bits are new.
t0 = time.time()
abl = [{"position_encoding": True, "seed": int(r["init"]), "test MAE": round(r["test MAE"], 2)}
for _, r in grid[grid["mask"] == 1].iterrows()]
for tm_seed in range(3):
m = fit_crtm(x_train, y_train, n_clauses=N_CLAUSES_W, T=15, weighted=True,
epochs=EPOCHS_W, batch_size=BATCH_W, seed=tm_seed, pos=False)
abl.append({"position_encoding": False, "seed": tm_seed,
"test MAE": round(mae(m, x_test, y_test), 2)})
print(f"({time.time() - t0:.0f}s)")
print(pd.DataFrame(abl).pivot(index="position_encoding", columns="seed",
values="test MAE").to_string())
(The True row is mask 1 of the section 7 grid — the same data and seeds — so only the three
runs without position bits are new.)
Without position literals every clause is a statement about some patch of the image, and since a 2×2 filter visits nine locations, "some patch has a 1 in its top-left corner" is true almost always. The clauses cannot isolate individual pixels and the error stays large — several hundred MAE, against 0 with position bits.
This is the same flag that
CLAUDE.mdwarns about from the other direction: on translation-invariant tasks (make_shapes,make_2d_noisy_xor) pinning clauses to locations hurts. Here the target is defined by absolute pixel positions, so it is essential. The flag is a statement about your data, not a tuning knob.
9. Noisy data, and the CNN comparison (Tables 6 and 7)¶
The other half of the paper's 72 datasets adds noise to the training data only, then
compares C-RTM against two CNNs: CNN1 with 50 filters and CNN2 with 1000. The paper does
not state the noise magnitude; the reported training MAEs of ≈5 for every model are consistent
with additive label noise of U(−10, 10) (mean absolute deviation 5), which is what we use
here. Test targets stay clean, which is why test MAEs come out lower than training MAEs — as
they do in the paper.
NOISE = 10.0
xs, y_clean, y_noisy, mask_n = make_weighted_mask_dataset(4000, SIZE, K, seed=1, noise=NOISE)
xtr_n, ytr_n, ytr_clean = xs[:3000], y_noisy[:3000], y_clean[:3000]
xte_n, yte_n = xs[3000:], y_clean[3000:]
print(f"label noise U(-{NOISE:g}, {NOISE:g}) on the training targets; "
f"mean |noise| = {(ytr_n - ytr_clean).abs().mean():.2f}")
class SmallCNN(nn.Module):
"""The paper's CNN1 / CNN2: a convolutional layer of n_filters followed by a linear head."""
def __init__(self, n_filters, size=SIZE, kernel=2):
super().__init__()
self.conv = nn.Conv2d(1, n_filters, kernel)
out = (size - kernel + 1) ** 2 * n_filters
self.head = nn.Linear(out, 1)
def forward(self, x):
return self.head(torch.relu(self.conv(x)).flatten(1)).squeeze(1)
def fit_cnn(n_filters, epochs=300, bs=64, lr=1e-3, seed=0):
torch.manual_seed(seed)
net = SmallCNN(n_filters).to(device)
opt = torch.optim.Adam(net.parameters(), lr=lr)
for _ in range(epochs):
perm = torch.randperm(xtr_n.shape[0], device=device)
for i in range(0, len(perm), bs):
j = perm[i:i + bs]
opt.zero_grad()
# train on targets scaled to [0, 1]: raw 0..1500 targets make the MSE badly
# conditioned and neither CNN converges at a sane learning rate.
nn.functional.mse_loss(net(xtr_n[j]), ytr_n[j] / Y_MAX).backward()
opt.step()
net.eval()
return net
EPOCHS_NOISY = 100
t0 = time.time()
rows = []
for name, nf in (("CNN1 (50 filters)", 50), ("CNN2 (1000 filters)", 1000)):
net = fit_cnn(nf)
with torch.no_grad():
rows.append({"model": name,
"train MAE (noisy targets)": round((net(xtr_n) * Y_MAX - ytr_n).abs().mean().item(), 2),
"test MAE (clean targets)": round((net(xte_n) * Y_MAX - yte_n).abs().mean().item(), 2)})
for label, n_clauses in (("C-RTM (15 weighted clauses)", 15),
("C-RTM (300 weighted clauses)", 300)):
m = fit_crtm(xtr_n, ytr_n, n_clauses=n_clauses, T=n_clauses, weighted=True,
epochs=EPOCHS_NOISY, seed=0)
if n_clauses == 15:
crtm_noisy = m
rows.append({"model": label,
"train MAE (noisy targets)": round(mae(m, xtr_n, ytr_n), 2),
"test MAE (clean targets)": round(mae(m, xte_n, yte_n), 2)})
print(f"({time.time() - t0:.0f}s)\n")
print(pd.DataFrame(rows).to_string(index=False))
print(f"\nnoise floor — the training MAE of a model that fits the clean function exactly: "
f"{(ytr_n - ytr_clean).abs().mean():.2f}")
The CNNs behave as the paper describes: both land near the noise floor, and CNN1 tests better than CNN2 despite having twenty times fewer filters — the paper's observation that "testing MAEs of CNN1 are lower than the MAEs of CNN2 […] a possible reason could be the overfitting of data by CNN2".
The C-RTM row does not reproduce the paper. Tables 6 and 7 report C-RTM at ≈5 training MAE
and ≈1.2 test MAE, in line with the CNNs. torchtsetlin gets nowhere near that (roughly 150 MAE
with 15 clauses here), and enlarging the clause budget — the paper's own prescription for noisy
data ("the C-RTM can be initialized with a large number of clauses (= T) […] the error
decreases") — narrows the gap (roughly 60 MAE with 300 clauses) without closing it. These runs
keep the short schedule of earlier versions of this notebook (100 epochs, batch 8), not section
4's longer one.
A plausible reading, consistent with section 7 but not separately tested: noise-free, the
machine has to land on one exact clause configuration. With noisy labels that configuration is
never quite right — every example still triggers feedback with probability
((y − ŷ)/ŷmax)² — so even a model that has found the clean function keeps being nudged, and
each 100-unit output step is a whole clause-weight unit; there is no "slightly wrong" to settle
into.
Things that were tried, with the older feedback rules of torchtsetlin 0.1.x, and did not
close the gap: sequential feedback and batch sizes of 1, 2 and 4 (so this is not
torchtsetlin's batched-feedback approximation); clause budgets of 60, 150, 300 and 1000; and
weighting by duplication instead of integer weights. The paper's T and s for these runs are
not stated — it inherits them from the Regression TM paper — so a hyper-parameter gap is the
most likely remaining explanation. Treat the noisy-data column as an open reproduction
question, not a settled result.
What survives is the asymmetry the paper closes on: the C-RTM's clauses are the printable table from section 5, and CNN2's thousand filters are not.
print("What the C-RTM learned from noisy labels (15 clauses; the clauses that fire):\n")
print(clause_conditions(crtm_noisy, xte_n).query("fires > 0").to_string(index=False))
print("\nimportant pixels and their place values:",
{tuple(p): int(mask_n[tuple(p)]) for p in (mask_n > 0).nonzero().tolist()})
print("\ncontribution recovered by probing the model (true values in brackets):")
rec_n = sensitivity_map(crtm_noisy, xte_n)
for p in (mask_n > 0).nonzero().tolist():
print(f" pixel {tuple(p)}: {rec_n[tuple(p)]:7.1f} [{100 * mask_n[tuple(p)]:.0f}]")
10. What this does — and does not — give you for segmentation¶
The C-RTM is in the segmentation folder because it is the only published TM mechanism for continuous, patch-derived output. Taking stock of what actually transfers:
Transfers.
- Clause outputs can be summed rather than argmax-ed, so a TM can emit a magnitude.
- Position literals give clauses absolute spatial addressing — the prerequisite for any output that is indexed by location.
- Integer clause weights make the "how many clauses vote" arithmetic compact and, as section 5 showed, directly readable.
Does not transfer.
- The output is one number per image, not per pixel. Equation 5 sums over all clauses and all patch positions; the location information is consumed to decide whether a clause fires, then discarded. Nothing in the C-RTM produces a map.
- The evaluation is entirely synthetic: 3×3 to 5×5 images, already-binary pixels (so booleanization is skipped), and 16 attainable output values.
- The exact-recovery result depends on knowing the right clause count and output resolution in advance — which you only know because the data was constructed — and, with the current library, on a clause count below the paper's and a long training schedule.
Turning this into a per-pixel output is the subject of the companion notebook,
05_ctm_unet_segmentation.ipynb, which works through the
one paper that actually attempts it.
Takeaways¶
ConvRegressionTsetlinMachineis the paper's model:y_range=(0, ŷmax)withT = 15reproduces Equation 5 exactly. The feedback probability is the square of Equation 6 withK = 1— the reference implementations' rule — which corrects small errors much more slowly.- Zero MAE is reproducible on noise-free data with integer clause weights — but on the current library only with fewer clauses than the paper's (8 rather than 15) and a longer schedule (400 epochs). With the paper's 15 clauses no restart reached it in our tests. The cause is the two rules adopted in 0.5 (weights floored at 1, squared feedback); with the 0.1.x rules the paper's settings worked. The paper's weighting-by-duplication did not converge here.
- The machine recovers the weight-mask, place value for place value — but you read it off by probing the model, not off the clause weights. The exact solution is not unique, and runs typically find case-split decompositions rather than the paper's canonical one.
- Exact recovery is a basin, not a limit — reached on about two thirds of the restarts here (8 of 12 in our run), never on the mask whose important pixels share one 2×2 window, and missed by whole 100-unit steps when missed. Selecting restarts by training MAE finds it without touching test labels whenever a restart found it (3 of 4 masks in our run).
position_encoding=Trueis load-bearing whenever the target depends on absolute location, and harmful when it doesn't.- The noisy-data result did not reproduce. The CNN comparison behaves as the paper
describes;
torchtsetlin's C-RTM lands far above the noise floor, and a larger clause budget narrows but does not close the gap. Section 9 says what was tried.