#pragma mode( separator(.,;) integer(h32) )

// ============================================================
// TRENDTEST - beta 3d
// Chi-square for linear trend in proportions (2 x k table)
// ============================================================
// Author: Roger Metcalf, 2026.  Written with AI assistance
//   (Anthropic Claude); concept, direction and all hardware
//   testing by the author.
//
// FREE FOR EDUCATIONAL AND PERSONAL USE.  Keep this header with
//   the file and credit the author if you pass it on.  Not for
//   sale or commercial redistribution.  Provided as is, with no
//   warranty.
// DISCLAIMER: a teaching and demonstration aid. Not validated for
//   forensic, clinical or legal use; do not rely on its results for
//   legal purposes.
//
// Firmware: written for the 2026-09-09 firmware (2.4.2).  This
//   version (beta 3d) has run on the PC harness only (ppl_run_v3e,
//   verify_TrendTest_beta3d.py): NOT on the emulator, NOT on
//   hardware.  No version has run on hardware yet.  Beta 2 ran on
//   the HP Prime Virtual Calculator in round 1 (17 Oct 2026) and
//   passed.  Beta 3b ran there in round 2 (17 Oct 2026): the CAS
//   answered page 4's calls, then trRepl stopped the program with
//   "Error: Invalid input"; the proportions plot read "(+9.6E-2
//   SE)"; Esc in Help turned the page.  Beta 3c fixed all three and
//   passed the retest that afternoon.  The call form (beta 3a) is
//   the one CASProbe beta 1b found on the HP Prime Virtual
//   Calculator (build 2026-09-09, 17 Oct 2026).
//
// beta 3d (19 Oct 2026) - what changed: the fleet review of 19 Oct
//   2026 (findings_C: improvement 1 as a fix, and C-T1 to C-T4).
//   Line numbers are beta 3c's.  Nothing else changed.
//   - TR_P KEEPS ITS DIGITS (trCompute, 762).  The trend p was
//     1 - CHISQUARE_CDF(1, X2).  The Prime's reals hold 12 digits, so
//     once the true tail is below about 5E-13 (X2 above about 50) the
//     CDF rounds to exactly 1 and 1 - CDF is exactly 0: on a 12-digit
//     model of the Prime the snoring table's TR_P is 0 (true p
//     1.5E-17); the calculator itself has not been asked yet (the
//     CHANGES' device checklist does).  For 1 df the
//     upper tail is a LOWER normal tail, which keeps its digits:
//       P(chi-square(1) > X2) = P(|Z| > sqrt(X2)) = 2 NORMALD_CDF(0, 1, -sqrt(X2)).
//     trPT now comes from that, so TR_P and every displayed trend p
//     (page 1, page 2's panel, the plot footer, the summary returned
//     to Home) are computed the same way.  The screens still print
//     "p < 0.0001" below 0.0001 (trPStr), as before.  Pearson's and
//     the departure's p (k-1 and k-2 df) have no such shortcut and
//     stay 1 - CDF.
//   - C-T1, the header's residual sentence (246-247).  The squared
//     residuals add up to Pearson - Armitage's N r^2; the departure
//     the screens show is Pearson - Mantel's (N-1) r^2, which is
//     larger by trend/(N-1).  The sentence now says so (no screen
//     changed; Mantel's form stays the headline).
//   - C-T2, help page 3 (trHelp, 1964-1968).  It listed R's
//     prop.trend.test and SAS PROC FREQ / TREND as "the same test
//     elsewhere" without saying they print Armitage's N r^2: 72.688
//     for the snoring table, where this program shows Mantel's
//     72.659.  The five lines now say which program prints which
//     form.  JMP is no longer named: which form it prints was not
//     checked.
//   - C-T3, the main menu (TrendTest, 2010).  "Plot: proportions +
//     trend line" (30 characters) was at the edge of the box that
//     cut CICatcher's 35-character line after about 31; it is now
//     "Plot: proportions, trend line" (29).
//   - C-T4, Enter data (trEnter, 876-892).  Cases or totals with a
//     fraction (24.5 of 100.5) were accepted: the table printed them
//     rounded ("25/101") beside a p computed from 24.5.  They are now
//     refused with a message: they count people.
//   - The DISCLAIMER line under the distribution terms (Roger, 19 Oct
//     2026).  No page of the program states the terms, so no screen
//     changed.
//
// beta 3c (17 Oct 2026) - what changed: the second emulator round's
//   findings (HP Prime Virtual Calculator 2.4.2, build 2026-09-09).
//   Line numbers are beta 3b's.  Nothing else changed.
//   - PAGE 4 NO LONGER STOPS (trRepl, 571-584).  On results page 4
//     the emulator said "Error: Invalid input" and the program
//     stopped.  INSTRING("", "a") is that error on the Prime (typed in
//     Home).  trRepl looked for the pattern again in what was left of
//     the text, and solve's answer {[a,b]} ends in the "}" it removes,
//     so what was left was "".  It now looks only while text is left,
//     and an empty pattern gives the text back unchanged.
//   - NOR DOES THE CAS FALLBACK (trCas 533, trCxClean 1090-1094,
//     trCxSol 1109-1122).  When the CAS fails, trCas sets its answer
//     to "" and then tested INSTRING(answer, "undef"): the same error.
//     The tests now run only on a non-empty answer (the kit's shape),
//     trCxClean hands "" back as "", and trCxSol stops when nothing is
//     left once the brackets are gone, or after "a=".  The other
//     INSTRINGs never see an empty text: trCxMake (1148) has returned
//     if trRat's text is empty, trRatNote (566) gets it only with a
//     "/", trCxSol's tests on the halves follow a length check, and
//     trCxFrac (1214) searches the literal "0123456789".
//   - NO SCIENTIFIC NUMBERS.  On the Prime STRING(x, 2, d) prints
//     scientific form when 0 < |x| < 10^-d (STRING(0.096, 2, 1) is
//     "9.6E-2", the E a small capital); the proportions plot read
//     "line 0.017 (+9.6E-2 SE)".  A new helper, trFx(x, d), rounds
//     first, so the value is 0 or at least 10^-d and the form stays
//     fixed.  Every display that can get such a value uses it: the
//     plot's readout and statistics (1453-1464), page 1's fit, p, CI,
//     the three statistics and the WLS line (915-930), page 2's
//     statistics (1008, 1027, 1031), the decimals beside the exact
//     fractions and on the decimal page (trDec 1224-1228, trNum
//     1233-1242) and the summary returned to Home (1928).  A score
//     below 0.01 now shows three decimals (trXStr 381-388), and an
//     odds ratio below 0.001 reads "< 0.001" (trOrStr 374-379), as
//     trPStr does for p.  At the other end, a sum of 1E12 or more on
//     the decimal page (midranks on N = 20000) is scientific on the
//     Prime whatever the decimals: trNum now keeps six significant
//     digits there (7.00018E13, not 7.0E13).
//   - ESC LEAVES HELP (trWaitAny 430-444, trHelp 1819-1861).  Esc
//     turned the help page like any other key.  trWaitAny now returns
//     the key, Esc goes back to the menu from any of the three pages,
//     and the footers say so: "any key = page 2, Esc = menu".
//
// beta 3b (17 Oct 2026) - what changed: the first emulator round's
//   findings (Emulator_Round1_Results_2026-10-17_v1.md).  Line numbers
//   are beta 3a's.  Nothing else changed.
//   - SHORTER EXAMPLES LINES (trExamples, 815-818).  CHOOSE cuts a line
//     wider than its box: the emulator showed "(made up" with the ")"
//     gone.  "Sleep hours: a U shape (made up)" is now "Sleep: a U
//     shape (made up)", and "Convictions, re-arrest (made up)" is
//     "Prior convictions (made up)" (the data's own label, shown with
//     the results, still says "convictions, re-arrest").
//   - SOLID DOTS (trDot, 379-386).  GfxProbe beta 1 (T09) showed that
//     ARC_P circles of radius 3, 2 and 1 plus the centre pixel leave
//     white holes on the diagonals: every dot shows a white X.  The
//     dots on the proportions and odds-ratio plots are now three
//     overlapping RECT_P bars, exactly a disc of radius 3, in the same
//     colours (one drawing call fewer per dot).  The black ring round
//     the selected level (two ARC_P circles, r 6 and 7) and the hollow
//     reference circle (r 3) are outlines on purpose and stay ARC_P.
//
// beta 3a (17 Oct 2026) - what changed: the CAS call form only.  CASProbe beta 1b
//   found that CAS(s) with s a variable hands the text back
//   UNEVALUATED, so beta 3's trCas (471-492) returned its own input;
//   the solve's echoed text kept "tra" in it, so trCxSol refused it,
//   and the page fell back to the decimals with the "decimal" badge.
//   Its IFERR line is now
//     STRING(EXPR("CAS(" + CHAR(34) + "string(" + s + ")" + CHAR(34) + ")"))
//   (BRIEF_CAS_v2), the quotes STRING() adds are dropped, and an
//   answer holding the infinity sign is refused like undef.  No
//   reader changed: trCxSol already drops braces, which is how the
//   Prime prints the solve of a system ({[23/10,4/5]}, CASProbe beta
//   1b T36), and the SSE and its derivatives are polynomials, which
//   print alike; its comment says so.  Nothing else changed.
//
// beta 3 (16 Oct 2026) - what changed: least squares as calculus,
//   built on beta 2; every beta 2 change is kept.  Line numbers
//   refer to beta 2.
//   - A FOURTH RESULTS PAGE, "Least squares as calculus" (Right
//     from page 3).  The line on the proportions plot is the
//     weighted least-squares line: it makes
//        SSE(a, b) = sum over the levels of n (p - a - b x)^2
//     (p = cases/n, x = the score in use) as small as it can be.
//     The page has the CAS write that SSE out for the table on the
//     screen, put it in order as a polynomial in a and b, take its
//     two partial derivatives, set them to 0 (the normal equations)
//     and solve them: the slope b and intercept a appear as exact
//     fractions beside the program's own decimals (trB, trA0).
//     The same line is ordinary least squares on all N people's
//     0/1 outcomes, and the trend chi-square asks whether b = 0.
//   - Only whole numbers go into CAS text: every count and score is
//     written as a whole number, or a whole number over a power of
//     ten (a midrank 1698.5 is 16985/10; the page says so).  A value
//     that needs more than six decimals cannot be written exactly,
//     and then the page shows the normal equations with the
//     program's decimal sums only.
//   - One call routine, trCas (the house CAS kit, with trIStr
//     writing whole numbers), so CASProbe's verdict is a one-line
//     change (made in beta 3a).  The CAS is asked when page 4 is shown and the counts
//     or scores have changed since it last answered -- never during
//     a drag.  A badge says "exact (CAS)", or "decimal" when the
//     CAS did not answer; the page then shows the decimal version.
//     trCas drops one pair of quotes round an answer (HP: answers
//     in numerical mode come back as strings); an answer with a
//     decimal point (approximate mode) counts as no answer.
//   - The results footers count four pages; help pages 2 and 3
//     describe page 4.  The kit's parameter names are trIv and trCs
//     (trN is TRN, the transpose, and trS is already the trend
//     sum S).  The drag (trPlotTouch) is unchanged.
// beta 2 (16 Oct 2026) - what changed: the touch and graphics
//   version, built on beta 1b; every beta 1b fix is kept.  Line
//   numbers refer to beta 1b.
//   - DRAG A DOT.  On the proportions plot a press within 10 px
//     of a level's dot, then an up-or-down drag, sets that level's
//     cases from the finger's height: cases = ROUND(p * n), held
//     between 0 and n (and so that the table keeps at least one
//     case and one non-case).  Every sample that changes the count
//     reruns trCompute and redraws: the line tilts, the Wilson bar
//     resizes, the footer's trend, Pearson and departure numbers
//     move.  The y-scale is frozen while the finger is down and
//     refitted on lift.  Up/Down do the same from the keyboard.
//     The data as entered are kept in trOX / trOA; the Restore
//     soft key (or 0) puts them back.
//   - SOFT KEYS {Prop, OR, Scores, Data, Restore, Menu} on both
//     plots and on the results pages (one view loop, trView, in
//     place of trResults 416-486 and trPlotLoop 669-706).  Scores
//     (or +) cycles the scores: as entered / 1,2,3... / midranks,
//     and recomputes -- the trend statistic moves, Pearson's does
//     not (help page 1's warning, 726-728, made live).  Data
//     reopens the entry with the current counts filled in.
//     Footer text now sits above y 218, under the soft keys'
//     strip; the plot box is 28..152 (was 30..186).
//   - A 95% confidence band round the WLS line (light magenta,
//     one pixel column at a time): fit +/- 1.96 SE with
//     SE^2 = P(1-P) (1/N + (x - xbar)^2 / Sxx), the weighted
//     regression of the proportions on the scores with weights
//     n / (P(1-P)) -- the variance the trend statistic uses.
//     Checked against statsmodels (verify_TrendTest.py).
//   - Dots coloured by their residual from the line: blue above,
//     red below, grey within 1 SE (SE^2 = P(1-P)/n); the selected
//     level is ringed in black (red now means "below the line").
//   - The WLS line is one clipped LINE_P (2 px) instead of 520
//     PIXON_P per frame (541-549), so a drag frame stays cheap.
//   - Entering the proportions plot animates the fit: the line
//     starts flat at the overall proportion and swings to the
//     fitted slope over 15 frames while the trend statistic counts
//     up (any key or tap skips).
//   - Results page 1 gains a "fit" column: the cases the fitted
//     line predicts at each level (they add up to A).
//   - A new results page 2: the chi-square(1) curve with the trend
//     statistic and its tail shaded, the chi-square(k-2) curve
//     with the departure statistic, and "1 df concentrates the
//     evidence" with the two 5% cut-offs.  Each density is cached
//     per pixel column.  The odds-ratio page is now page 3.
//   - Menu -> Examples (a CHOOSE): snoring as before, plus two
//     ILLUSTRATIVE data sets with made-up counts -- a U-shaped
//     sleep-hours example where the departure test fires, and a
//     criminal-justice one (prior convictions 0/1/2/3+ against
//     re-arrest within a year).  Made-up data are labelled in red
//     on the results page and in the plot titles.
//   - The OR plot's footer says why its axis is logarithmic.
//   - Help is three pages (soft keys, drag, keys, band, colours).
//   - trCompute no longer resets the highlighted level (348), so a
//     drag keeps its dot; it also sets the x-window, Sxx and
//     P(1-P).  Scores print with 0, 1 or 2 decimals each (trXStr)
//     instead of one setting for all (trXDec, 338-342).
// beta 1b (16 Oct 2026) - what changed: fixes from the fleet review
//   of 16 Oct; no change to the arithmetic.  Line numbers refer to
//   beta 1a.
//   - Two names the Prime reads as built-in commands (it ignores
//     case): the totals list trN was TRN (transpose) and doT in
//     trDrawOR was DOT -- the same fault that stopped ChiTable's
//     first paste (csC = CSC).  Now trTot and doTk.
//   - The plot-top constant trPT0 is now trPTop (98): it was one
//     typo away from trPT, the p-value (92).
//   - The menu items carry no numbers (CHOOSE numbers them itself)
//     (742-747); the "no data yet" message (755) and help page 2
//     (701, 704) name the items instead of quoting numbers.
//   - Footers in the house form (541, 625, 627; 422, 455, 694,
//     713).
//   - MAX(trHi) on a list (480) is replaced by a running maximum
//     taken inside the Wilson loop (279-287), kept in trHiMx.
//   - A tap selects a level only inside the plot box (656-662).
// beta 1a (11 Oct 2026): clean exits.  Every "wait for the finger
//   to lift" loop now goes through trDrain(), which gives up after
//   one second, so a phantom touch cannot freeze the program; and
//   Quit returns the summary line to Home instead of PRINTing it,
//   so no Terminal screen is left open.  The arithmetic was
//   checked line for line against scipy and statsmodels (300
//   random tables) before the PPL was written.
//
// Run TrendTest() from the Program Catalog.
// ============================================================
//
// THE QUESTION
//   k exposure levels have a natural ORDER and a score x_i
//   (a dose, hours per day, 0/1/2/3 ...).  At each level a_i of
//   n_i people have the outcome.  Pearson's chi-square (k-1 df)
//   only asks "do the proportions differ?".  The chi-square for
//   trend (1 df) asks "do they rise or fall STEADILY with the
//   score?" -- the question an epidemiologist usually means.
//
// THE ARITHMETIC (Mantel 1963; Epi Info StatCalc's form)
//   N = sum n_i,  A = sum a_i  (all cases),  P = A / N
//   S = sum(x_i a_i) - A * sum(x_i n_i) / N
//   V = A (N-A) / (N^2 (N-1)) * [ N sum(x_i^2 n_i) - (sum x_i n_i)^2 ]
//   X2 trend = S^2 / V,  1 df   ( = (N-1) r^2, r = Pearson
//     correlation of score with the 0/1 outcome over all N )
//   Armitage (1955) uses N r^2 instead; the difference is
//   N/(N-1), invisible at classroom sample sizes.  R's
//   prop.trend.test and SAS PROC FREQ / TREND print Armitage's
//   N r^2 (72.688 for the snoring table, where this program shows
//   Mantel's 72.659; beta 3d, fleet review C-T2).
//   X2 Pearson for the 2 x k table, k-1 df.
//   X2 departure = X2 Pearson - X2 trend, k-2 df: large means the
//     pattern is not a straight line.
//   Weighted least-squares line through the proportions
//     (weights n_i): slope b = S / Sxx, Sxx = sum n_i (x_i - xbar)^2,
//     xbar = sum x_i n_i / N; the line passes through (xbar, P).
//   Its 95% band: fit(x) +/- 1.96 * sqrt( P(1-P) * (1/N +
//     (x - xbar)^2 / Sxx) ) -- the weighted regression with the
//     variance of each p_i taken as P(1-P)/n_i, the null variance
//     the trend statistic itself uses: (b / SE_b)^2 = N r^2.
//   Residual of a level: p_i - fit_i, in units of
//     SE_i = sqrt(P(1-P)/n_i); the squared residuals add up to
//     Pearson - Armitage's N r^2.  The departure statistic the
//     screens show is Pearson - Mantel's (N-1) r^2, which is larger
//     by trend/(N-1): 53.1149 against the residuals' 53.1141 for the
//     sleep example (beta 3d, fleet review C-T1).
//   Midrank scores: level i gets the average rank its n_i people
//     would have if all N were ranked by level (n_1 + ... +
//     n_(i-1) + (n_i + 1)/2).  With midranks the trend test is the
//     Wilcoxon-Mann-Whitney test in disguise.
//   Each level: p_i with a 95% Wilson score interval; the odds
//     ratio against level 1 with a Woolf (log) 95% interval.
//     If any cell is 0, 0.5 is added to every cell (Haldane-
//     Anscombe) and the screens say so.
//
// EXAMPLES (Menu -> Examples); every number below was checked
//   with SciPy and statsmodels (verify_TrendTest.py).
//   1 Snoring and heart disease, Norton & Dunn 1985, used by
//     Agresti.  Scores 0, 2, 4, 5; cases 24, 35, 21, 30 of 1379,
//     638, 213, 254.  X2 trend 72.659 (p < 0.0001), Pearson
//     72.782 (3 df), departure 0.123 (2 df, p = 0.940); slope
//     0.0200 per unit of score; ORs vs level 1: 3.277 (1.932-
//     5.558), 6.175 (3.373-11.306), 7.561 (4.341-13.172).
//     Scores 1,2,3,4: trend 72.008, departure 0.774 (p = 0.679);
//     midranks: trend 68.190, departure 4.592 (p = 0.101);
//     Pearson 72.782 every time.
//   2 Sleep hours and illness -- ILLUSTRATIVE, MADE-UP COUNTS (the
//     U shape echoes published sleep studies; the numbers do not
//     come from one).  Hours 5, 6, 7, 8, 9; cases 38, 40, 36, 44,
//     44 of 160, 320, 480, 320, 160.  X2 trend 1.105 (p = 0.293),
//     Pearson 54.220 (4 df, p < 0.0001), departure 53.115 (3 df,
//     p < 0.0001): no straight-line trend, but the levels differ --
//     high at both ends.  Residuals +4.1, -0.4, -4.1, -0.6, +4.3 SE.
//   3 Prior convictions and re-arrest within a year --
//     ILLUSTRATIVE, MADE-UP COUNTS.  Convictions 0, 1, 2, 3+ (3+
//     scored 3); re-arrested 72, 66, 51, 88 of 400, 240, 150, 200.
//     X2 trend 47.719 (p < 0.0001), Pearson 47.967 (3 df),
//     departure 0.248 (2 df, p = 0.883); slope 0.0852 per
//     conviction; ORs vs none: 1.728 (1.180-2.530), 2.347 (1.537-
//     3.584), 3.579 (2.452-5.225).  Midranks: trend 47.022.
//
// SCREENS
//   Enter data: how many levels (2-6), then the scores, the cases
//     and the totals in three dialogs.
//   Results, four pages (left/right arrows or a tap on the left/
//     right half): 1 the table with the fit column and the tests;
//     2 the two chi-square curves; 3 odds ratios and the formulas;
//     4 least squares as calculus (CAS; beta 3).
//   Plot: proportions with Wilson bars, the WLS trend line with
//     its 95% band, the overall proportion; drag a dot.
//   Plot: odds ratios vs level 1 on a log scale, Woolf bars.
//   Soft keys on all five: Prop, OR, Scores, Data, Restore, Menu.
//   Examples (a choice of three), Help (three pages).
//
// HOME EXPORTS (after a run)
//   TR_X2T  chi-square for trend     TR_P   its p-value (1 df;
//           beta 3d: a lower normal tail, so it keeps its digits
//           far out -- 1.5E-17 for snoring, not 0)
//   TR_X2   Pearson chi-square       TR_DF  its df (k - 1)
//   TR_OR   list of odds ratios vs level 1
//   TR_B    slope of the WLS line, per unit of score
// ============================================================


EXPORT TR_X2T, TR_P, TR_X2, TR_DF, TR_OR, TR_B;

// ---- the data ----
LOCAL trK, trX, trA, trTot, trHas, trLbl;
//  trK levels; trX scores in use, trA cases, trTot totals (lists;
//  trTot rather than trN, because trN is TRN, the transpose);
//  trHas = 1 once data are loaded; trLbl names the data set.
LOCAL trOX, trOA, trScM, trEdit, trLbl2, trIllu;
//  trOX, trOA: the scores and cases as entered (or as the example
//  gave them); Restore puts them back after a drag.
//  trScM: the scores in use -- 0 as entered, 1 equal steps 1,2,3...,
//  2 midranks.  trEdit = 1 once a drag or Up/Down changed a count.
//  trLbl2: a second label line; trIllu = 1 for made-up example data.

// ---- the results (set in trCompute) ----
LOCAL trNT, trAT, trS, trV, trX2T, trPT;
LOCAL trX2P, trPP, trX2D, trPD, trDFD, trB;
LOCAL trA0, trXbar, trP, trLo, trHi, trOR;
LOCAL trORL, trORH, trHald, trZ, trHiMx;
//  trHiMx: the largest Wilson upper limit, for the plot's y-scale
LOCAL trSxx, trPQ, trXlo, trXhi;
//  trSxx = sum n (x - xbar)^2 and trPQ = P (1 - P): the band and
//  the residual SEs need them.  trXlo, trXhi: the x-window both
//  plots and the touch handler use (the scores, padded by 10%).

// ---- plot box, shared by the drawing code and the touch handler ----
LOCAL trPL, trPR, trPTop, trPB, trSel;
//  trSel: the level currently highlighted on a plot
LOCAL trYhi, trYStp, trYDec, trFrz, trDX, trDY;
//  trYhi, trYStp, trYDec: plot 1's y-scale (top, tick step, tick
//  decimals); trFrz = 1 freezes it while a finger drags a dot.
//  trDX, trDY: each dot's pixel centre as last drawn -- the hit
//  boxes the touch handler reads, so drawing and touch agree.

// ---- house colours, and the chi-square glyph ----
LOCAL trBLUE, trMAG, trRED, trGRN, trORNG, trGRAY;
LOCAL trFOOT, trGRID, trLORN, trLBLU, trLMAG, trChi;

// ---- beta 3: least squares as calculus (results page 4) ----
LOCAL trCasOK, trCasMx, trCxK, trCxSt, trCxNt, trCxE;
LOCAL trCxA, trCxB, trCxSa, trCxSb;
//  trCasOK 1 = the CAS answered ("exact (CAS)" badge), 0 = it did not
//          ("decimal"); trCasMx the longest CAS text the program sends
//  trCxK   the SSE text the CAS last answered for (when the table or
//          the scores change, so does this text, and the CAS is asked
//          again); trCxSt what happened: 1 exact answers, 2 the CAS
//          did not answer, 3 a value cannot be written exactly, 4 too
//          long to send, 0 not asked yet
//  trCxNt  an example of a value written over a power of ten, for the
//          page ("1698.5 = 16985/10"), or ""
//  trCxE   the SSE as the CAS put it in order; trCxA, trCxB its two
//          partial derivatives; trCxSa, trCxSb the solution a and b


// ------------------------------------------------------------
// SMALL HELPERS
// ------------------------------------------------------------

// A "nice" step (1, 2 or 5 times a power of ten) at least as big
// as the value handed in -- Heckbert's rule, for axis ticks.
trNice(nnV)
BEGIN
  LOCAL nnE, nnF;
  IF nnV <= 0 THEN RETURN 1; END;
  nnE := FLOOR(LOG(nnV));
  nnF := nnV / 10^nnE;
  IF nnF <= 1 THEN RETURN 10^nnE; END;
  IF nnF <= 2 THEN RETURN 2 * 10^nnE; END;
  IF nnF <= 5 THEN RETURN 5 * 10^nnE; END;
  RETURN 10^(nnE + 1);
END;

// The sum of a list (a FOR loop, so nothing outside the verified
// command set is needed).
trSum(smL)
BEGIN
  LOCAL smI, smT;
  smT := 0;
  FOR smI FROM 1 TO SIZE(smL) DO
    smT := smT + smL(smI);
  END;
  RETURN smT;
END;

// A fresh copy of a list, built element by element.  On the Prime
// "L2 := L1" already copies, but writing the copy out keeps the
// saved originals (trOA) safe whatever the list semantics -- the
// PC test bench shares a list between two names after ":=".
trCopy(cyL)
BEGIN
  LOCAL cyI, cyR;
  cyR := {};
  FOR cyI FROM 1 TO SIZE(cyL) DO
    cyR := CONCAT(cyR, {cyL(cyI)});
  END;
  RETURN cyR;
END;

// beta 3c: fixed decimals that stay fixed.  On the Prime STRING(x, 2, d)
// prints SCIENTIFIC form when 0 < |x| < 10^-d: STRING(0.096, 2, 1) is
// "9.6E-2" with a small-capital E (U+1D07; Home, emulator, 17 Oct 2026),
// and the proportions plot showed "(+9.6E-2 SE)".  ROUND first makes the
// value either 0 or at least 10^-d, so the fixed form always holds.  A
// value that rounds to 0 prints 0.000, never -0.000.
trFx(fxV, fxD)
BEGIN
  LOCAL fxR;
  fxR := ROUND(fxV, fxD);
  IF fxR == 0 THEN fxR := 0; END;
  RETURN STRING(fxR, 2, fxD);
END;

// A p-value as text.  Tiny values are reported as "< 0.0001":
// four decimals are all a screen needs, and Pearson's and the
// departure's p, which are 1 - CDF, lose their digits far out anyway.
// (beta 3d: the trend's p is a lower normal tail and keeps its
// digits -- see TR_P in Home -- but it prints the same way here.)
trPStr(psP)
BEGIN
  IF psP < 0.0001 THEN RETURN "p < 0.0001"; END;
  RETURN "p = " + STRING(psP, 2, 4);
END;

// Decimals that show a coefficient's first three significant
// digits (at least 2, at most 6).
trDecOf(dcV)
BEGIN
  IF dcV == 0 THEN RETURN 2; END;
  RETURN MAX(2, MIN(6, 2 - FLOOR(LOG(ABS(dcV)))));
END;

// An odds ratio as text: 3 decimals below 10, 2 below 100, else 1.
// beta 3c: below 0.001 it says "< 0.001", as trPStr does for p.  Three
// decimals printed it in scientific form on the Prime (a level dragged
// to no cases against a level of all cases gives OR 0.0000003).
trOrStr(osV)
BEGIN
  IF osV >= 100 THEN RETURN STRING(osV, 2, 1); END;
  IF osV >= 10 THEN RETURN STRING(osV, 2, 2); END;
  IF osV < 0.001 THEN RETURN "< 0.001"; END;
  RETURN STRING(osV, 2, 3);
END;

// A score as text: no decimals for a whole number, one for a
// half (midranks are whole or end in .5), two otherwise.
// beta 3c: three, rounded first (trFx), for a score below 0.01 (a dose
// of 0.005): two decimals printed it in scientific form on the Prime.
trXStr(xsV)
BEGIN
  IF FP(xsV) == 0 THEN RETURN STRING(xsV, 2, 0); END;
  IF FP(2 * xsV) == 0 THEN RETURN STRING(xsV, 2, 1); END;
  IF ABS(xsV) >= 0.01 THEN RETURN STRING(xsV, 2, 2); END;
  RETURN trFx(xsV, 3);
END;

// The name of the scores in use (trScM).
trScName()
BEGIN
  IF trScM == 1 THEN RETURN "1,2,3..."; END;
  IF trScM == 2 THEN RETURN "midranks"; END;
  RETURN "as entered";
END;

// A filled dot of radius 3 (7 px across).  beta 3b: three overlapping
// RECT_P bars -- 7 x 3, 3 x 7 and 5 x 5 -- which together cover
// exactly the pixels of a disc of radius 3.  Up to beta 3a this was
// three concentric ARC_P circles and the centre pixel, which leave
// holes on the diagonals (GfxProbe T09, emulator 17 Oct 2026).
trDot(dtG, dtX, dtY, dtCol)
BEGIN
  RECT_P(dtG, dtX - 3, dtY - 1, dtX + 3, dtY + 1, dtCol, dtCol);
  RECT_P(dtG, dtX - 1, dtY - 3, dtX + 1, dtY + 3, dtCol, dtCol);
  RECT_P(dtG, dtX - 2, dtY - 2, dtX + 2, dtY + 2, dtCol, dtCol);
END;

// A dotted horizontal line (LINE_P has no dash style).
trDotLine(dlG, dlX1, dlX2, dlY, dlCol)
BEGIN
  LOCAL dlI;
  FOR dlI FROM dlX1 TO dlX2 STEP 3 DO
    PIXON_P(dlG, dlI, dlY, dlCol);
  END;
END;

// Wait for the finger to lift -- but never for more than a second.
// MOUSE(1) is -1 while nothing touches the screen.  Without the
// time limit a stuck or phantom touch would hold the program here
// for ever, and the only way out would be a reset.
trDrain()
BEGIN
  LOCAL drT;
  drT := TICKS;                                  // milliseconds
  REPEAT UNTIL MOUSE(1) == -1 OR TICKS - drT > 1000;
END;

// Wait for any key or tap (drain, wait, drain).  The help pages
// use it; the results pages and plots have their own loop (trView).
// beta 3c: returns the key code (-1 for a tap), so help can leave on
// Esc (4) from any page.
trWaitAny()
BEGIN
  LOCAL waK, waM, waM1, waR;
  REPEAT waK := GETKEY; UNTIL waK == -1;
  trDrain();
  REPEAT
    waK := GETKEY;
    waM := MOUSE;
    waM1 := waM(1);
  UNTIL waK > -1 OR SIZE(waM1) > 0;
  waR := waK;
  trDrain();
  REPEAT waK := GETKEY; UNTIL waK == -1;
  RETURN waR;
END;

// Show the frame composed in G1, then the soft keys.  DRAWMENU
// paints the SCREEN, not G1, so it comes after the blit and is
// repeated every frame (ChiTable, G2 and G1).  A tap with y >= 218
// is a soft key: slot MIN(6, IP(x / 53.4) + 1).  The labels never
// produce a key code, so the hardware keys stay as shortcuts.
trShow()
BEGIN
  BLIT_P(G0, 0, 0, G1);
  DRAWMENU({"Prop", "OR", "Scores", "Data", "Restore", "Menu"});
END;

// Ask for akM numbers (2 to 6) in one INPUT dialog.  INPUT needs
// variable NAMES, so the values go through the locals a1..a6 and
// come back as a list.  akD is the list of defaults.  Returns {}
// if the user cancels.
trAsk(akM, akTtl, akLbl, akHlp, akD)
BEGIN
  LOCAL a1, a2, a3, a4, a5, a6;
  LOCAL akOK, akP;
  akP := CONCAT(akD, {0, 0, 0, 0, 0, 0});
  a1 := akP(1); a2 := akP(2); a3 := akP(3); a4 := akP(4); a5 := akP(5); a6 := akP(6);
  akOK := 0;
  CASE
    IF akM == 2 THEN akOK := INPUT({a1, a2}, akTtl, {akLbl + " 1:", akLbl + " 2:"}, {akHlp, akHlp}); END;
    IF akM == 3 THEN akOK := INPUT({a1, a2, a3}, akTtl, {akLbl + " 1:", akLbl + " 2:", akLbl + " 3:"}, {akHlp, akHlp, akHlp}); END;
    IF akM == 4 THEN akOK := INPUT({a1, a2, a3, a4}, akTtl, {akLbl + " 1:", akLbl + " 2:", akLbl + " 3:", akLbl + " 4:"}, {akHlp, akHlp, akHlp, akHlp}); END;
    IF akM == 5 THEN akOK := INPUT({a1, a2, a3, a4, a5}, akTtl, {akLbl + " 1:", akLbl + " 2:", akLbl + " 3:", akLbl + " 4:", akLbl + " 5:"}, {akHlp, akHlp, akHlp, akHlp, akHlp}); END;
    IF akM == 6 THEN akOK := INPUT({a1, a2, a3, a4, a5, a6}, akTtl, {akLbl + " 1:", akLbl + " 2:", akLbl + " 3:", akLbl + " 4:", akLbl + " 5:", akLbl + " 6:"}, {akHlp, akHlp, akHlp, akHlp, akHlp, akHlp}); END;
  END;
  IF akOK == 0 THEN RETURN {}; END;
  RETURN SUB({a1, a2, a3, a4, a5, a6}, 1, akM);
END;


// ------------------------------------------------------------
// THE CAS KIT (the house shapes, BRIEF_CAS 16 Oct 2026; beta 3)
// beta 3a: trCas in the call form CASProbe beta 1b found on the
// emulator (BRIEF_CAS_v2), not yet run on a calculator.  The kit's
// parameters are
// trIv and trCs here: "trN" would be TRN (the transpose, case is
// ignored) and trS is already the trend sum S.
// ------------------------------------------------------------

// A whole number as plain digits: no digit-grouping marks, no decimal point, no
// exponent.  Home's STRING() follows the Home number format (PyProbe saw grouping
// marks and a special minus sign), which the CAS cannot read, so CAS text is built
// only from this routine.
trIStr(trIv)
BEGIN
  LOCAL dS, dV, dD;
  dV := ROUND(ABS(trIv), 0);
  IF dV == 0 THEN RETURN "0"; END;
  dS := "";
  WHILE dV > 0 DO
    dD := dV MOD 10;
    dS := CHAR(48 + dD) + dS;
    dV := (dV - dD) / 10;
  END;
  IF trIv < 0 THEN dS := "-" + dS; END;
  RETURN dS;
END;

// The ONE place this program calls the CAS.  It returns the CAS answer as text,
// or "" if the CAS failed, gave undef, an error or an infinity.
// The call form is the one CASProbe beta 1b found on the HP Prime Virtual
// Calculator (build 2026-09-09, 17 Oct 2026): CAS(s) with s a variable hands
// the text back UNEVALUATED, so the text is built into a literal CAS("...")
// that EXPR reads at run time, and the CAS is asked for string(...) of the
// answer, so a big integer crosses into PPL whole, not as a 12-digit real.
trCas(trCs)
BEGIN
  LOCAL cR;
  cR := "";
  IFERR
    cR := STRING(EXPR("CAS(" + CHAR(34) + "string(" + trCs + ")" + CHAR(34) + ")"));
  THEN
    cR := "";
  END;
  // The answer is a PPL string, and STRING() of it adds one pair of quotes
  // (CASProbe beta 1b T12, emulator): drop that pair.
  IF DIM(cR) >= 2 THEN
    IF LEFT(cR, 1) == CHAR(34) AND MID(cR, DIM(cR), 1) == CHAR(34) THEN
      cR := MID(cR, 2, DIM(cR) - 2);
    END;
  END;
  // undef, an error, or an infinity (CHAR(8734) is the infinity sign; it
  // also catches the Prime's plus-minus infinity, 1/0): no usable answer
  // an empty answer is already a failure; INSTRING("", ...) is "Invalid input" on
  // the Prime (emulator, 17 Oct 2026), so the checks run only on a non-empty answer
  IF DIM(cR) > 0 THEN
    IF INSTRING(cR, "undef") > 0 OR INSTRING(cR, "rror") > 0 OR INSTRING(cR, CHAR(8734)) > 0 THEN cR := ""; END;
  END;
  RETURN cR;
END;

// A number as EXACT CAS text: a whole number ("24", "(-1)"), or a
// whole number over a power of ten ("(16985/10)" for 1698.5) -- or
// "" when it needs more than six decimals (1/3, say) or is too big
// to write exactly.  Built with trIStr only, so no decimal point and
// no Home number format ever reaches the CAS.
trRat(rtV)
BEGIN
  LOCAL rtD, rtS, rtN;
  FOR rtD FROM 0 TO 6 DO
    rtS := rtV * 10^rtD;
    rtN := ROUND(rtS, 0);
    IF ABS(rtS - rtN) <= 0.000001 + ABS(rtS) * 0.000000000001 THEN
      IF ABS(rtN) >= 100000000000 THEN RETURN ""; END;
      IF rtD == 0 THEN
        IF rtN < 0 THEN RETURN "(" + trIStr(rtN) + ")"; END;
        RETURN trIStr(rtN);
      END;
      RETURN "(" + trIStr(rtN) + "/" + trIStr(10^rtD) + ")";
    END;
  END;
  RETURN "";
END;

// For the page: "1698.5 = 16985/10" -- a value and the exact text
// trRat wrote for it (rnR holds a "/"; its denominator's zeros give
// the decimals to show).
trRatNote(rnV, rnR)
BEGIN
  LOCAL rnP, rnD;
  rnP := INSTRING(rnR, "/");
  rnD := MID(rnR, rnP + 1, DIM(rnR) - rnP - 1);
  RETURN STRING(rnV, 2, DIM(rnD) - 1) + " = " + MID(rnR, 2, DIM(rnR) - 2);
END;

// Every rpO in the text rpT replaced by rpN (INSTRING, LEFT, MID).
// beta 3c: INSTRING("", s) is "Error: Invalid input" on the Prime
// (emulator, 17 Oct 2026), and a pattern at the very end -- the "}"
// of solve's {[a,b]} -- leaves rpS empty, which stopped page 4.  So
// the position is looked for only while text is left; an empty
// pattern gives the text back unchanged (INSTRING(s, "") is 1, which
// would loop for ever).
trRepl(rpT, rpO, rpN)
BEGIN
  LOCAL rpR, rpS, rpI, rpK;
  IF DIM(rpO) == 0 THEN RETURN rpT; END;
  rpR := ""; rpS := rpT; rpK := DIM(rpO);
  rpI := 0; IF DIM(rpS) > 0 THEN rpI := INSTRING(rpS, rpO); END;
  WHILE rpI > 0 DO
    IF rpI > 1 THEN rpR := rpR + LEFT(rpS, rpI - 1); END;
    rpR := rpR + rpN;
    IF rpI + rpK > DIM(rpS) THEN rpS := ""; ELSE rpS := MID(rpS, rpI + rpK); END;
    rpI := 0; IF DIM(rpS) > 0 THEN rpI := INSTRING(rpS, rpO); END;
  END;
  RETURN rpR + rpS;
END;

// Cut a long formula into lines of at most wpW characters, at the
// best place in the second half of each line: just after a closing
// bracket or a times sign, just before a + or - that joins two
// terms, or just after a / that is followed by a bracket -- and of
// those, the one with the fewest brackets open.  A number or a
// fraction such as 63080529409/7933174634 is never split; only a
// formula with no such place is cut at the width.  Returns the list
// of lines.
trWrap(wpT, wpW)
BEGIN
  LOCAL wpL, wpS, wpI, wpB, wpC, wpD;
  LOCAL wpBd, wpP, wpN, wpOk;
  wpL := {}; wpS := wpT;
  WHILE DIM(wpS) > wpW DO
    wpB := 0; wpBd := 999; wpD := 0;
    FOR wpI FROM 1 TO wpW DO
      wpC := MID(wpS, wpI, 1);
      wpP := ""; IF wpI > 1 THEN wpP := MID(wpS, wpI - 1, 1); END;
      wpN := MID(wpS, wpI + 1, 1);
      // a break BEFORE a + or - (one that joins terms, not a sign)
      IF wpI > IP(wpW / 2) AND (wpC == "+" OR wpC == "-") THEN
        IF wpP <> "(" AND wpP <> "^" AND wpP <> CHAR(215) AND wpP <> "/" AND wpD <= wpBd THEN
          wpB := wpI - 1; wpBd := wpD;
        END;
      END;
      IF wpC == "(" THEN wpD := wpD + 1; END;
      IF wpC == ")" THEN wpD := wpD - 1; END;
      // a break AFTER a ), a times sign, or a / followed by a bracket
      wpOk := 0;
      IF wpC == ")" OR wpC == CHAR(215) THEN wpOk := 1; END;
      IF wpC == "/" AND wpN == "(" THEN wpOk := 1; END;
      IF wpI > IP(wpW / 2) AND wpOk == 1 AND wpD <= wpBd THEN
        wpB := wpI; wpBd := wpD;
      END;
    END;
    IF wpB == 0 THEN wpB := wpW; END;
    wpL := CONCAT(wpL, {LEFT(wpS, wpB)});
    wpS := MID(wpS, wpB + 1);
  END;
  IF DIM(wpS) > 0 THEN wpL := CONCAT(wpL, {wpS}); END;
  RETURN wpL;
END;


// ------------------------------------------------------------
// THE ARITHMETIC
// ------------------------------------------------------------

// The scores in use, from trScM: the scores as entered (trOX), equal
// steps 1, 2, 3 ..., or midranks.  The trend statistic depends on
// the SPACING of the scores; Pearson's chi-square does not use them.
trSetX()
BEGIN
  LOCAL sxI, sxC;
  IF trScM == 0 THEN trX := trCopy(trOX); END;
  IF trScM == 1 THEN
    trX := {};
    FOR sxI FROM 1 TO trK DO
      trX := CONCAT(trX, {sxI});
    END;
  END;
  IF trScM == 2 THEN
    // midrank: the average rank the n_i people at level i would get
    // if all N were lined up by level (ties share their ranks)
    trX := {};
    sxC := 0;                                    // people at the levels below
    FOR sxI FROM 1 TO trK DO
      trX := CONCAT(trX, {sxC + (trTot(sxI) + 1) / 2});
      sxC := sxC + trTot(sxI);
    END;
  END;
END;

trCompute()
BEGIN
  LOCAL cpI, cpSxa, cpSxn, cpSxxn, cpEa, cpEb;
  LOCAL cpP, cpDen, cpCen, cpHalf, cpSxx, cpAdj;
  LOCAL cpA1, cpB1, cpAi, cpBi, cpO, cpSe;
  LOCAL cpPad;

  // totals and the three sums the trend statistic needs
  trNT := 0; trAT := 0; cpSxa := 0; cpSxn := 0; cpSxxn := 0;
  FOR cpI FROM 1 TO trK DO
    trNT := trNT + trTot(cpI);
    trAT := trAT + trA(cpI);
    cpSxa := cpSxa + trX(cpI) * trA(cpI);
    cpSxn := cpSxn + trX(cpI) * trTot(cpI);
    cpSxxn := cpSxxn + trX(cpI)^2 * trTot(cpI);
  END;

  // chi-square for trend (Mantel): S^2 / V, 1 df.  V > 0 because
  // trEnter, the drag and Up/Down all keep 0 < A < N.
  trS := cpSxa - trAT * cpSxn / trNT;
  trV := trAT * (trNT - trAT) / (trNT^2 * (trNT - 1)) * (trNT * cpSxxn - cpSxn^2);
  trX2T := trS^2 / trV;
  // beta 3d: the trend p as a LOWER tail, so it keeps its digits.
  // With 1 df, chi-square = Z^2, so P(chi-square > X2) is
  // P(|Z| > sqrt(X2)) = 2 * P(Z < -sqrt(X2)).  The Prime's reals hold
  // 12 digits: once the true tail is below about 5E-13 (X2 above about
  // 50), CHISQUARE_CDF rounds to exactly 1 and the old
  // 1 - CHISQUARE_CDF(1, X2) is exactly 0 -- on a 12-digit model of the
  // Prime the snoring table's TR_P is 0 instead of 1.5E-17.  A lower
  // tail is computed directly and loses nothing (fleet review 19 Oct
  // 2026, TrendTest improvement 1).  trX2T >= 0 (a square over V > 0),
  // and ^0.5 is the house square root.
  trPT := 2 * NORMALD_CDF(0, 1, -(trX2T^0.5));

  // Pearson chi-square for the 2 x k table, k-1 df
  trX2P := 0;
  FOR cpI FROM 1 TO trK DO
    cpEa := trTot(cpI) * trAT / trNT;            // expected cases
    cpEb := trTot(cpI) - cpEa;                   // expected non-cases
    trX2P := trX2P + (trA(cpI) - cpEa)^2 / cpEa + (trTot(cpI) - trA(cpI) - cpEb)^2 / cpEb;
  END;
  trPP := 1 - CHISQUARE_CDF(trK - 1, trX2P);

  // departure from linearity, k-2 df (needs at least 3 levels)
  trDFD := trK - 2;
  trX2D := MAX(0, trX2P - trX2T);
  IF trDFD >= 1 THEN
    trPD := 1 - CHISQUARE_CDF(trDFD, trX2D);
  ELSE
    trPD := -1;
  END;

  // weighted least-squares line through the proportions
  trXbar := cpSxn / trNT;
  cpSxx := cpSxxn - trNT * trXbar^2;
  trB := trS / cpSxx;
  trA0 := trAT / trNT - trB * trXbar;
  // what the band and the residual colours need: Sxx and P(1-P)
  trSxx := cpSxx;
  trPQ := (trAT / trNT) * (1 - trAT / trNT);
  // the x-window of both plots: the scores, padded by 10% each side
  cpPad := (trX(trK) - trX(1)) * 0.1;
  trXlo := trX(1) - cpPad;
  trXhi := trX(trK) + cpPad;

  // each level: proportion and its 95% Wilson score interval
  trP := {}; trLo := {}; trHi := {};
  trHiMx := 0;                                   // running maximum of trHi
  FOR cpI FROM 1 TO trK DO
    cpP := trA(cpI) / trTot(cpI);
    cpDen := 1 + trZ^2 / trTot(cpI);
    cpCen := (cpP + trZ^2 / (2 * trTot(cpI))) / cpDen;
    cpHalf := trZ * (cpP * (1 - cpP) / trTot(cpI) + trZ^2 / (4 * trTot(cpI)^2))^0.5 / cpDen;
    trP := CONCAT(trP, {cpP});
    trLo := CONCAT(trLo, {MAX(0, cpCen - cpHalf)});
    trHi := CONCAT(trHi, {MIN(1, cpCen + cpHalf)});
    IF trHi(cpI) > trHiMx THEN trHiMx := trHi(cpI); END;
  END;

  // odds ratios against level 1, with Woolf 95% intervals.
  // A zero cell would make an OR 0 or infinite, so then 0.5 is
  // added to every cell (Haldane-Anscombe) and trHald says so.
  trHald := 0;
  FOR cpI FROM 1 TO trK DO
    IF trA(cpI) == 0 OR trA(cpI) == trTot(cpI) THEN trHald := 1; END;
  END;
  cpAdj := 0.5 * trHald;
  cpA1 := trA(1) + cpAdj;
  cpB1 := trTot(1) - trA(1) + cpAdj;
  trOR := {}; trORL := {}; trORH := {};
  FOR cpI FROM 1 TO trK DO
    cpAi := trA(cpI) + cpAdj;
    cpBi := trTot(cpI) - trA(cpI) + cpAdj;
    cpO := (cpAi * cpB1) / (cpBi * cpA1);
    cpSe := (1 / cpAi + 1 / cpBi + 1 / cpA1 + 1 / cpB1)^0.5;
    trOR := CONCAT(trOR, {cpO});
    trORL := CONCAT(trORL, {EXP(LN(cpO) - trZ * cpSe)});
    trORH := CONCAT(trORH, {EXP(LN(cpO) + trZ * cpSe)});
  END;

  // Home exports.  (beta 1b also reset the highlighted level here;
  // beta 2 leaves it, so a dragged dot stays selected -- trEnter
  // and the examples set it instead.)
  TR_X2T := trX2T; TR_P := trPT; TR_X2 := trX2P; TR_DF := trK - 1;
  TR_OR := trOR; TR_B := trB;
  trHas := 1;
END;


// ------------------------------------------------------------
// DATA ENTRY
// ------------------------------------------------------------
trEnter()
BEGIN
  LOCAL enK, enX, enA, enN, enI, enBad;
  LOCAL enDX, enDA, enDN;
  enK := 0;
  CHOOSE(enK, "How many exposure levels?",
    {"2 levels", "3 levels", "4 levels", "5 levels", "6 levels"});
  IF enK == 0 THEN RETURN 0; END;
  enK := enK + 1;

  // the dialogs open on the data already loaded when it has this
  // many levels (so the Data soft key edits what is on the screen);
  // otherwise on scores 1, 2, 3 ... and zero counts
  IF trHas == 1 AND enK == trK THEN
    enDX := trOX; enDA := trA; enDN := trTot;
  ELSE
    enDX := MAKELIST(K, K, 1, enK, 1);
    enDA := MAKELIST(0, K, 1, enK, 1);
    enDN := MAKELIST(0, K, 1, enK, 1);
  END;

  // scores: 1, 2, 3 ... unless real doses are known
  enX := trAsk(enK, "Scores, in order of exposure", "Score",
    "A number per level: 0,1,2.. or the dose itself", enDX);
  IF SIZE(enX) == 0 THEN RETURN 0; END;

  enA := trAsk(enK, "Cases (with the outcome)", "Cases",
    "How many at this level have the outcome", enDA);
  IF SIZE(enA) == 0 THEN RETURN 0; END;

  enN := trAsk(enK, "Totals (everyone at the level)", "Total n",
    "How many people are at this level in all", enDN);
  IF SIZE(enN) == 0 THEN RETURN 0; END;

  // ---- checks that teach ----
  // beta 3d: cases and totals count people, so a fraction (24.5 cases
  // of 100.5) is refused (enBad 4).  Before, the arithmetic ran on
  // 24.5 while the table printed "25/101" (fleet review 19 Oct 2026,
  // C-T4).  Scores may have fractions: a dose of 0.5 is fine.
  enBad := 0;
  FOR enI FROM 1 TO enK DO
    IF enN(enI) <= 0 THEN enBad := 1; END;
    IF enA(enI) < 0 OR enA(enI) > enN(enI) THEN enBad := 2; END;
    IF enI > 1 THEN
      IF enX(enI) <= enX(enI - 1) THEN enBad := 3; END;
    END;
    IF FP(enA(enI)) <> 0 OR FP(enN(enI)) <> 0 THEN enBad := 4; END;
  END;
  IF enBad == 1 THEN MSGBOX("Every total n must be above 0."); RETURN 0; END;
  IF enBad == 2 THEN MSGBOX("Cases must be between 0 and the total at that level."); RETURN 0; END;
  IF enBad == 3 THEN MSGBOX("Scores must increase from level 1 to level " + STRING(enK, 2, 0) + ": the levels are ordered, and the order is the point of the test."); RETURN 0; END;
  IF enBad == 4 THEN MSGBOX("Cases and totals count people: whole numbers only."); RETURN 0; END;
  IF trSum(enA) == 0 OR trSum(enA) == trSum(enN) THEN
    MSGBOX("Every level has the same outcome (all cases, or none), so there is no trend to test.");
    RETURN 0;
  END;
  IF trSum(enN) < 2 THEN MSGBOX("At least two people are needed."); RETURN 0; END;

  trK := enK; trTot := enN;
  trOX := trCopy(enX);                           // kept for Restore
  trOA := trCopy(enA);
  trA := trCopy(enA);
  trScM := 0; trEdit := 0; trIllu := 0; trSel := 1;
  trSetX();
  trLbl := "Your data: " + STRING(trK, 2, 0) + " levels, N = " + STRING(trSum(enN), 2, 0) + ", " + STRING(trSum(enA), 2, 0) + " cases";
  trLbl2 := "";
  trCompute();
  RETURN 1;
END;

// The example data sets.  Returns 1 if one was loaded, 0 on Cancel.
//   1 Snoring and heart disease (Norton & Dunn 1985, in Agresti):
//     real data, a clean straight-line trend.
//   2 Sleep hours and illness: MADE-UP counts, U-shaped, so the
//     trend test sees almost nothing and the departure test fires.
//   3 Prior convictions and re-arrest within a year: MADE-UP
//     counts; "3+" is scored 3 -- try the midranks with Scores.
trExamples()
BEGIN
  LOCAL exI;
  exI := 0;
  // beta 3b: lines shortened -- CHOOSE cuts a line wider than its box
  // (the emulator showed "(made up" cut); every line is now 27 or fewer
  CHOOSE(exI, "Examples", {
    "Snoring and heart disease",
    "Sleep: a U shape (made up)",
    "Prior convictions (made up)"});
  IF exI == 0 THEN RETURN 0; END;
  CASE
    IF exI == 1 THEN
      trK := 4;
      trOX := {0, 2, 4, 5};
      trOA := {24, 35, 21, 30};
      trTot := {1379, 638, 213, 254};
      trLbl := "Example: snoring and heart disease";
      trLbl2 := "never, occasional, nearly every night, every night";
      trIllu := 0;
    END;
    IF exI == 2 THEN
      trK := 5;
      trOX := {5, 6, 7, 8, 9};
      trOA := {38, 40, 36, 44, 44};
      trTot := {160, 320, 480, 320, 160};
      trLbl := "ILLUSTRATIVE (made-up counts): sleep, illness";
      trLbl2 := "hours of sleep a night: 5, 6, 7, 8, 9";
      trIllu := 1;
    END;
    DEFAULT
      trK := 4;
      trOX := {0, 1, 2, 3};
      trOA := {72, 66, 51, 88};
      trTot := {400, 240, 150, 200};
      trLbl := "ILLUSTRATIVE (made up): convictions, re-arrest";
      trLbl2 := "prior convictions 0, 1, 2, 3+ (3+ scored 3)";
      trIllu := 1;
  END;
  trA := trCopy(trOA);
  trScM := 0; trEdit := 0; trSel := 1;
  trSetX();
  trCompute();
  RETURN 1;
END;


// ------------------------------------------------------------
// RESULTS: three drawn pages (composed in G1; trShow puts them up
// with the soft keys).  The frames stop at y 217, above the strip.
// ------------------------------------------------------------

// Page 1: the table (with the cases the line predicts) and the tests
trPage1()
BEGIN
  LOCAL rsI, rsY, rsS, rsDec, rsF;
  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 217, trORNG, trLORN);
  TEXTOUT_P("Chi-square test for linear trend", G1, 10, 6, 2, trBLUE);
  IF trIllu == 1 THEN
    TEXTOUT_P(trLbl, G1, 10, 22, 1, trRED, 300);
  ELSE
    TEXTOUT_P(trLbl, G1, 10, 22, 1, trGRAY, 300);
  END;
  IF trEdit == 1 THEN
    TEXTOUT_P("Counts changed by dragging; Restore = original", G1, 10, 33, 1, trRED, 300);
  ELSE
    TEXTOUT_P(trLbl2, G1, 10, 33, 1, trGRAY, 300);
  END;
  TEXTOUT_P("level", G1, 10, 46, 1, trGRAY);
  TEXTOUT_P("score", G1, 42, 46, 1, trGRAY);
  TEXTOUT_P("cases/n", G1, 80, 46, 1, trGRAY);
  TEXTOUT_P("fit", G1, 134, 46, 1, trMAG);
  TEXTOUT_P("p", G1, 176, 46, 1, trGRAY);
  TEXTOUT_P("95% CI", G1, 218, 46, 1, trGRAY);
  LINE_P(G1, 10, 57, 309, 57, trGRAY);
  FOR rsI FROM 1 TO trK DO
    rsY := 60 + (rsI - 1) * 12;
    // fit: the cases the weighted line predicts, n * (a + b x)
    rsF := trTot(rsI) * (trA0 + trB * trX(rsI));
    TEXTOUT_P(STRING(rsI, 2, 0), G1, 16, rsY, 1, #000000h);
    TEXTOUT_P(trXStr(trX(rsI)), G1, 42, rsY, 1, #000000h, 36);
    TEXTOUT_P(STRING(trA(rsI), 2, 0) + "/" + STRING(trTot(rsI), 2, 0), G1, 80, rsY, 1, #000000h, 52);
    // beta 3c: trFx, so a small fit, p or CI end stays in fixed decimals
    TEXTOUT_P(trFx(rsF, 1), G1, 134, rsY, 1, trMAG, 40);
    TEXTOUT_P(trFx(trP(rsI), 4), G1, 176, rsY, 1, trBLUE);
    TEXTOUT_P("(" + trFx(trLo(rsI), 3) + "-" + trFx(trHi(rsI), 3) + ")", G1, 218, rsY, 1, trBLUE);
  END;
  // the tests
  TEXTOUT_P("fit = cases the WLS line predicts (they sum to A)", G1, 10, 134, 1, trMAG, 300);
  TEXTOUT_P("Trend (1 df):  " + trChi + " = " + trFx(trX2T, 3) + "   " + trPStr(trPT), G1, 10, 146, 1, trMAG, 300);
  TEXTOUT_P("Pearson (" + STRING(trK - 1, 2, 0) + " df):  " + trChi + " = " + trFx(trX2P, 3) + "   " + trPStr(trPP), G1, 10, 158, 1, trBLUE, 300);
  IF trDFD >= 1 THEN
    TEXTOUT_P("Departure (" + STRING(trDFD, 2, 0) + " df):  " + trChi + " = " + trFx(trX2D, 3) + "   " + trPStr(trPD), G1, 10, 170, 1, trGRN, 300);
  ELSE
    TEXTOUT_P("Departure from linearity: needs 3+ levels.", G1, 10, 170, 1, trGRN, 300);
  END;
  rsDec := trDecOf(trB);
  IF trB < 0 THEN rsS := " - " + trFx(-trB, rsDec); ELSE rsS := " + " + trFx(trB, rsDec); END;
  TEXTOUT_P("WLS line:  p = " + trFx(trA0, rsDec) + rsS + " * score", G1, 10, 182, 1, trGRAY, 300);
  TEXTOUT_P("Scores: " + trScName() + "  (the Scores key cycles them)", G1, 10, 194, 1, trGRAY, 300);
  TEXTOUT_P("1/4, ⇆ or tap = page, Esc = menu", G1, 10, 205, 1, trFOOT, 300);
END;

// One small chi-square panel for page 2: the curve for cqDf df
// between pixel columns cqX1 and cqX2, the tail beyond the statistic
// cqSt shaded (its area is the p-value cqPv), the 5% cut-off dotted
// red, the statistic as a solid line in cqCol.
trChiPanel(cqX1, cqX2, cqDf, cqSt, cqPv, cqCol, cqNm)
BEGIN
  LOCAL cqW, cqDen, cqI, cqXmx, cqCr, cqYmx;
  LOCAL cqStp, cqDec, cqV, cqN, cqPx, cqPy;
  LOCAL cqLast, cqSx, cqCx, cqT, cqB;
  cqT := 34; cqB := 120;                         // the panel's top and bottom
  cqW := cqX2 - cqX1;
  cqCr := CHISQUARE_ICDF(cqDf, 0.95);            // the 5% cut-off
  // x from 0 to past the cut-off, far enough to include the
  // statistic -- unless it is more than 4 cut-offs out; then it
  // is marked at the right edge as off the scale
  cqXmx := MAX(cqCr * 1.6, MIN(cqSt * 1.15, cqCr * 4));
  cqStp := trNice(cqXmx / 4);
  cqXmx := CEILING(cqXmx / cqStp - 0.000001) * cqStp;
  cqDec := MAX(0, -FLOOR(LOG(cqStp)));
  // the height: df 1 rises without limit near 0, so it is simply
  // clipped; df 2 starts at 0.5; larger df peak at df - 2
  cqYmx := 0.6;
  IF cqDf == 2 THEN cqYmx := 0.55; END;
  IF cqDf > 2 THEN cqYmx := CHISQUARE(cqDf, cqDf - 2) * 1.15; END;
  // ---- the density, once per pixel column (Lessons Learned 6):
  // the shaded tail and the curve both read cqDen ----
  cqDen := {};
  FOR cqI FROM 1 TO cqW DO
    cqDen := CONCAT(cqDen, {CHISQUARE(cqDf, cqI / cqW * cqXmx)});
  END;
  TEXTOUT_P(cqNm + ": " + trChi + "(" + STRING(cqDf, 2, 0) + ")", G1, cqX1, 22, 1, cqCol, cqW);
  // ---- x ticks and grid ----
  cqV := 0; cqN := 0;
  WHILE cqV <= cqXmx + cqStp / 1000 AND cqN < 10 DO
    cqPx := cqX1 + IP(cqV / cqXmx * cqW + 0.5);
    LINE_P(G1, cqPx, cqT, cqPx, cqB, trGRID);
    LINE_P(G1, cqPx, cqB, cqPx, cqB + 3, #000000h);
    TEXTOUT_P(STRING(cqV, 2, cqDec), G1, cqPx - 5, cqB + 5, 1, trGRAY, 30);
    cqV := cqV + cqStp;
    cqN := cqN + 1;
  END;
  // ---- the p-value: the tail beyond the statistic, column by column ----
  cqSx := cqX1 + IP(MIN(cqSt, cqXmx) / cqXmx * cqW + 0.5);
  IF cqSt < cqXmx THEN
    FOR cqPx FROM MAX(cqSx, cqX1 + 1) TO cqX2 DO
      cqPy := cqB - IP(MIN(cqDen(cqPx - cqX1) / cqYmx, 1) * (cqB - cqT));
      IF cqPy < cqB THEN LINE_P(G1, cqPx, cqB - 1, cqPx, cqPy, trLMAG); END;
    END;
  END;
  // ---- the curve, pixel column by pixel column ----
  cqLast := -1;
  FOR cqPx FROM cqX1 + 1 TO cqX2 DO
    cqPy := cqB - IP(MIN(cqDen(cqPx - cqX1) / cqYmx, 1) * (cqB - cqT));
    IF cqLast >= 0 THEN LINE_P(G1, cqPx - 1, cqLast, cqPx, cqPy, trBLUE); END;
    cqLast := cqPy;
  END;
  LINE_P(G1, cqX1, cqT, cqX1, cqB, #000000h);
  LINE_P(G1, cqX1, cqB, cqX2, cqB, #000000h);
  // ---- the 5% cut-off, dotted red ----
  cqCx := cqX1 + IP(cqCr / cqXmx * cqW + 0.5);
  FOR cqI FROM cqT TO cqB STEP 3 DO
    PIXON_P(G1, cqCx, cqI, trRED);
  END;
  TEXTOUT_P("5%", G1, cqCx + 3, cqT + 14, 1, trRED, 24);
  // ---- the statistic: a solid line, or at the edge if off the scale ----
  IF cqSt <= cqXmx THEN
    LINE_P(G1, cqSx, cqT, cqSx, cqB, cqCol);
    LINE_P(G1, cqSx + 1, cqT, cqSx + 1, cqB, cqCol);
  ELSE
    LINE_P(G1, cqX2, cqT, cqX2, cqB, cqCol);
    LINE_P(G1, cqX2 - 1, cqT, cqX2 - 1, cqB, cqCol);
    TEXTOUT_P("further " + CHAR(8594), G1, cqX2 - 52, cqT + 2, 1, cqCol, 50);
  END;
  TEXTOUT_P(trChi + " = " + trFx(cqSt, 2) + ", " + trPStr(cqPv), G1, cqX1, 140, 1, cqCol, cqW + 6);
END;

// Page 2: where the two statistics fall
trPageChi()
BEGIN
  LOCAL pcCt, pcCp;
  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  TEXTOUT_P("Where the two statistics fall", G1, 10, 4, 2, trBLUE);
  trChiPanel(14, 154, 1, trX2T, trPT, trMAG, "trend");
  IF trDFD >= 1 THEN
    trChiPanel(170, 310, trDFD, trX2D, trPD, trGRN, "departure");
  ELSE
    TEXTOUT_P("departure: needs", G1, 176, 60, 1, trGRN, 130);
    TEXTOUT_P("3 or more levels", G1, 176, 72, 1, trGRN, 130);
  END;
  pcCt := CHISQUARE_ICDF(1, 0.95);
  pcCp := CHISQUARE_ICDF(trK - 1, 0.95);
  IF trDFD >= 1 THEN
    TEXTOUT_P("Pearson " + trChi + " " + trFx(trX2P, 2) + " = trend " + trFx(trX2T, 2) + " + departure " + trFx(trX2D, 2), G1, 10, 158, 1, #000000h, 300);
    TEXTOUT_P("1 df concentrates the evidence: its 5% cut-off", G1, 10, 172, 1, trMAG, 300);
    TEXTOUT_P("is " + STRING(pcCt, 2, 2) + ", against " + STRING(pcCp, 2, 2) + " for Pearson's " + STRING(trK - 1, 2, 0) + " df.", G1, 10, 184, 1, trMAG, 300);
  ELSE
    TEXTOUT_P("Pearson " + trChi + " = " + trFx(trX2P, 2) + " on 1 df: with 2 levels", G1, 10, 158, 1, #000000h, 300);
    TEXTOUT_P("Pearson and trend ask the same question", G1, 10, 172, 1, trMAG, 300);
    TEXTOUT_P("(they differ only by N/(N-1)).", G1, 10, 184, 1, trMAG, 300);
  END;
  TEXTOUT_P("2/4, ⇆ or tap = page, Esc = menu", G1, 10, 205, 1, trFOOT, 300);
END;

// Page 3: odds ratios and the formulas
trPage3()
BEGIN
  LOCAL rsI, rsY;
  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 217, trMAG, trLMAG);
  TEXTOUT_P("Odds ratios against level 1", G1, 10, 6, 2, trMAG);
  TEXTOUT_P("level", G1, 10, 26, 1, trGRAY);
  TEXTOUT_P("score", G1, 54, 26, 1, trGRAY);
  TEXTOUT_P("OR", G1, 104, 26, 1, trGRAY);
  TEXTOUT_P("95% CI (Woolf)", G1, 164, 26, 1, trGRAY);
  LINE_P(G1, 10, 37, 309, 37, trGRAY);
  FOR rsI FROM 1 TO trK DO
    rsY := 40 + (rsI - 1) * 12;
    TEXTOUT_P(STRING(rsI, 2, 0), G1, 16, rsY, 1, #000000h);
    TEXTOUT_P(trXStr(trX(rsI)), G1, 54, rsY, 1, #000000h, 46);
    IF rsI == 1 THEN
      TEXTOUT_P("1", G1, 104, rsY, 1, trMAG);
      TEXTOUT_P("(reference level)", G1, 164, rsY, 1, trGRAY);
    ELSE
      TEXTOUT_P(trOrStr(trOR(rsI)), G1, 104, rsY, 1, trMAG);
      TEXTOUT_P("(" + trOrStr(trORL(rsI)) + " - " + trOrStr(trORH(rsI)) + ")", G1, 164, rsY, 1, trMAG);
    END;
  END;
  IF trHald == 1 THEN
    TEXTOUT_P("Zero cell: +0.5 in every cell (Haldane).", G1, 10, 114, 1, trRED, 300);
  END;
  TEXTOUT_P(trChi + " trend = S^2 / V  (1 df),  where", G1, 10, 128, 1, #000000h);
  TEXTOUT_P("S = sum(x a) - A sum(x n) / N", G1, 10, 140, 1, #000000h);
  TEXTOUT_P("V = A (N-A) / (N^2 (N-1))", G1, 10, 152, 1, #000000h);
  TEXTOUT_P("      * [ N sum(x^2 n) - (sum x n)^2 ]", G1, 10, 164, 1, #000000h);
  TEXTOUT_P("Departure = Pearson - trend, k-2 df", G1, 10, 176, 1, #000000h);
  TEXTOUT_P("Home: TR_X2T TR_P TR_X2 TR_DF TR_OR TR_B", G1, 10, 190, 1, trGRAY, 300);
  TEXTOUT_P("3/4, ⇆ or tap = page, Esc = menu", G1, 10, 205, 1, trFOOT, 300);
END;


// ------------------------------------------------------------
// PAGE 4 (beta 3): LEAST SQUARES AS CALCULUS
// The line on the proportions plot is the weighted least-squares
// line: of all lines p = a + b x it makes
//     SSE(a, b) = sum over the levels of  n (p - a - b x)^2
// smallest (p = cases/n, x = the score in use, weights n).  At the
// bottom of that bowl both partial derivatives are 0 -- the two
// "normal equations" -- and solving them gives the program's
//     b = S / Sxx  and  a = P - b xbar.
// The CAS does the algebra on the table on screen, exactly.
// ------------------------------------------------------------

// A CAS answer that is no use as an EXACT answer: one with a
// decimal point (the CAS is in its approximate mode) counts as no
// answer, and the page falls back to the decimals.
// beta 3c: "" (the CAS did not answer) stays "" without INSTRING, which
// stops the program on an empty text ("Invalid input", emulator).
trCxClean(ccR)
BEGIN
  IF DIM(ccR) == 0 THEN RETURN ""; END;
  IF INSTRING(ccR, ".") > 0 THEN RETURN ""; END;
  RETURN ccR;
END;

// Read a and b out of solve's answer, whatever brackets the CAS
// puts round them -- {[a,b]} on the Prime (CASProbe beta 1b T36,
// emulator 17 Oct 2026), [[a,b]] or list[[a,b]] in desktop Giac, or
// [a=...,b=...] -- into trCxSa and trCxSb; both stay "" if the
// answer is not exactly two values.
trCxSol(soR)
BEGIN
  LOCAL soS, soK;
  trCxSa := ""; trCxSb := "";
  IF DIM(soR) == 0 THEN RETURN; END;
  soS := trRepl(soR, "list", "");
  soS := trRepl(soS, "[", ""); soS := trRepl(soS, "]", "");
  soS := trRepl(soS, "{", ""); soS := trRepl(soS, "}", "");
  soS := trRepl(soS, " ", "");
  // beta 3c: nothing left (solve gave {} or []): no answer.  INSTRING
  // on an empty text stops the program on the Prime (emulator).
  IF DIM(soS) == 0 THEN RETURN; END;
  soK := INSTRING(soS, ",");
  IF soK == 0 THEN soK := INSTRING(soS, ";"); END;
  IF soK < 2 OR soK >= DIM(soS) THEN RETURN; END;
  trCxSa := LEFT(soS, soK - 1);
  trCxSb := MID(soS, soK + 1);
  IF INSTRING(trCxSb, ",") > 0 OR INSTRING(trCxSb, ";") > 0 THEN
    trCxSa := ""; trCxSb := "";                  // more than two values
    RETURN;
  END;
  IF INSTRING(trCxSa, "=") > 0 THEN trCxSa := MID(trCxSa, INSTRING(trCxSa, "=") + 1); END;
  IF INSTRING(trCxSb, "=") > 0 THEN trCxSb := MID(trCxSb, INSTRING(trCxSb, "=") + 1); END;
  // beta 3c: "a=" with nothing after it leaves an empty text: no answer
  IF DIM(trCxSa) == 0 OR DIM(trCxSb) == 0 THEN trCxSa := ""; trCxSb := ""; RETURN; END;
  // a symbol left in it means no single solution: not usable
  IF INSTRING(trCxSa, "tr") > 0 OR INSTRING(trCxSb, "tr") > 0 THEN trCxSa := ""; trCxSb := ""; END;
END;

// Ask the CAS for page 4, if the table or the scores changed since
// it last answered.  The SSE is written out for the data on screen,
//   n1*(a1/n1-tra-trb*x1)^2 + n2*(a2/n2-tra-trb*x2)^2 + ...
// in the symbols tra (a) and trb (b), every number exact (trRat).
// Then four calls, through trCas only:
//   normal(SSE)                     the SSE in order, a polynomial
//   normal(diff(SSE,tra))           its slope along a
//   normal(diff(SSE,trb))           its slope along b
//   solve([d/da = 0, d/db = 0], [tra, trb])   the normal equations
trCxMake()
BEGIN
  LOCAL cxE, cxI, cxN, cxA, cxX, cxR;
  cxE := ""; trCxNt := "";
  FOR cxI FROM 1 TO trK DO
    cxN := trRat(trTot(cxI));
    cxA := trRat(trA(cxI));
    cxX := trRat(trX(cxI));
    IF DIM(cxN) == 0 OR DIM(cxA) == 0 OR DIM(cxX) == 0 THEN
      // a value with more than six decimals: no exact text
      trCxSt := 3; trCxK := ""; trCasOK := 0;
      RETURN;
    END;
    // one example of a value written over a power of ten, for the page
    IF DIM(trCxNt) == 0 AND INSTRING(cxX, "/") > 0 THEN
      trCxNt := trRatNote(trX(cxI), cxX);
    END;
    IF cxI > 1 THEN cxE := cxE + "+"; END;
    cxE := cxE + cxN + "*(" + cxA + "/" + cxN + "-tra-trb*" + cxX + ")^2";
  END;
  IF cxE == trCxK THEN RETURN; END;              // answered already
  trCxK := cxE;
  trCxE := ""; trCxA := ""; trCxB := ""; trCxSa := ""; trCxSb := "";
  trCasOK := 0;
  // the longest call (solve) holds the SSE twice
  IF 2 * DIM(cxE) + 60 > trCasMx THEN trCxSt := 4; RETURN; END;
  // the G1's CAS speed is not known: say what is happening
  RECT_P(G0, 60, 100, 260, 124, trFOOT, #FFFFFFh);
  TEXTOUT_P("asking the CAS ...", G0, 106, 107, 1, trFOOT);
  trCxE := trCxClean(trCas("normal(" + cxE + ")"));
  trCxA := trCxClean(trCas("normal(diff(" + cxE + ",tra))"));
  trCxB := trCxClean(trCas("normal(diff(" + cxE + ",trb))"));
  cxR := trCxClean(trCas("solve([diff(" + cxE + ",tra)=0,diff(" + cxE + ",trb)=0],[tra,trb])"));
  trCxSol(cxR);
  trCxSt := 2;
  IF DIM(trCxE) > 0 AND DIM(trCxA) > 0 AND DIM(trCxB) > 0 AND DIM(trCxSa) > 0 THEN
    trCxSt := 1; trCasOK := 1;
  END;
END;

// A CAS answer written for people: tra -> a, trb -> b, * -> times.
trCxShow(csT)
BEGIN
  LOCAL csS;
  csS := trRepl(csT, "tra", "a");
  csS := trRepl(csS, "trb", "b");
  csS := trRepl(csS, "*", CHAR(215));
  RETURN csS;
END;

// Draw clT wrapped at 49 characters from y clY, in colour clC, in at
// most clMx lines; if it needs more, say how long it is instead (an
// exact answer is never cut silently).  Returns the next free y.
trCxLines(clT, clY, clMx, clC, clNm)
BEGIN
  LOCAL clL, clI, clYy;
  clL := trWrap(clT, 49);
  clYy := clY;
  IF SIZE(clL) > clMx THEN
    TEXTOUT_P(clNm + ": " + STRING(DIM(clT), 2, 0) + " characters, too long to show", G1, 10, clYy, 1, trGRAY, 300);
    RETURN clYy + 11;
  END;
  FOR clI FROM 1 TO SIZE(clL) DO
    TEXTOUT_P(clL(clI), G1, 10, clYy, 1, clC, 302);
    clYy := clYy + 11;
  END;
  RETURN clYy;
END;

// One line: the exact fraction beside the program's decimal,
// "b = 95339/4757909   (program 0.020038)" -- or, when the fraction
// is too long for the line, how many digits it has (an exact answer
// is never cut silently).  Returns the next free y.
trCxFrac(cfNm, cfF, cfD, cfY)
BEGIN
  LOCAL cfI, cfN, cfS;
  cfS := cfNm + " = " + cfF + "   (program " + trDec(cfD) + ")";
  IF DIM(cfS) > 49 THEN
    cfN := 0;
    FOR cfI FROM 1 TO DIM(cfF) DO
      IF INSTRING("0123456789", MID(cfF, cfI, 1)) > 0 THEN cfN := cfN + 1; END;
    END;
    cfS := cfNm + ": exact fraction, " + STRING(cfN, 2, 0) + " digits (program " + trDec(cfD) + ")";
  END;
  TEXTOUT_P(cfS, G1, 10, cfY, 1, trMAG, 302);
  RETURN cfY + 11;
END;

// The program's own decimal beside an exact fraction: about seven
// significant digits (0.02003799, 0.0000528181, 16.63624).
// beta 3c: rounded first (trFx), so a value below 1E-10 -- a slope that
// is 0 but for the last digit -- prints "0", never in scientific form.
trDec(deV)
BEGIN
  LOCAL deD;
  IF deV == 0 THEN RETURN "0"; END;
  deD := MAX(2, MIN(10, 6 - FLOOR(LOG(ABS(deV)))));
  IF ROUND(deV, deD) == 0 THEN RETURN "0"; END;
  RETURN trFx(deV, deD);
END;

// A number for the decimal page: whole numbers as they are, others
// with about six significant digits and no trailing zeros (260.8,
// not 260.800).  MID, not RIGHT: RIGHT is not in the verified list.
// beta 3c: rounded first (trFx).  A sum below 1E-6 (scores such as -2/3,
// 1/3, 1/3 leave -1E-11 on the Prime) printed in scientific form, and
// the zero strip below then cut the exponent's last 0 (1E-10 -> 1E-1).
trNum(nmV)
BEGIN
  LOCAL nmS, nmD;
  IF FP(nmV) == 0 AND ABS(nmV) < 1000000000 THEN RETURN STRING(nmV, 2, 0); END;
  IF nmV == 0 THEN RETURN "0"; END;
  // beta 3c: from 1E12 up the Prime prints scientific form whatever the
  // decimals (it holds 12 digits; midranks on N = 20000 give sum n x^2 =
  // 7E13).  Five decimals keep six significant digits, and the zero
  // strip below must not reach the exponent.
  IF ABS(nmV) >= 1000000000000 THEN RETURN trFx(nmV, 5); END;
  nmD := MAX(1, MIN(6, 5 - FLOOR(LOG(ABS(nmV)))));
  IF ROUND(nmV, nmD) == 0 THEN RETURN "0"; END;
  nmS := trFx(nmV, nmD);
  WHILE DIM(nmS) > 1 AND MID(nmS, DIM(nmS), 1) == "0" DO nmS := LEFT(nmS, DIM(nmS) - 1); END;
  IF MID(nmS, DIM(nmS), 1) == "." THEN nmS := LEFT(nmS, DIM(nmS) - 1); END;
  RETURN nmS;
END;

// Page 4: least squares as calculus.  Exact (the CAS's lines) when
// trCxSt = 1; otherwise the normal equations with the program's own
// decimal sums.  Either way, the program's decimals at the bottom.
trPage4()
BEGIN
  LOCAL p4Y, p4I, p4Sx, p4Sxx, p4Sxa, p4S;
  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 217, trBLUE, trLBLU);
  TEXTOUT_P("Least squares as calculus", G1, 10, 6, 2, trBLUE);
  IF trCasOK == 1 AND trCxSt == 1 THEN
    TEXTOUT_P("exact (CAS)", G1, 244, 8, 1, trFOOT);
  ELSE
    TEXTOUT_P("decimal", G1, 268, 8, 1, trFOOT);
  END;
  TEXTOUT_P("The plotted line makes this SSE smallest:", G1, 10, 24, 1, #000000h, 302);
  TEXTOUT_P("SSE = sum n(p - a - bx)^2; p = cases/n, x = score", G1, 10, 35, 1, #000000h, 302);
  p4Y := 46;
  IF trCxSt == 1 THEN
    // ---- exact: the CAS's own lines ----
    p4Y := trCxLines("SSE = " + trCxShow(trCxE), p4Y, 3, trMAG, "SSE");
    TEXTOUT_P("Smallest where both slopes are 0 (CAS diff):", G1, 10, p4Y, 1, #000000h, 302);
    p4Y := p4Y + 11;
    p4Y := trCxLines("dSSE/da = " + trCxShow(trCxA) + " = 0", p4Y, 2, trBLUE, "dSSE/da = 0");
    p4Y := trCxLines("dSSE/db = " + trCxShow(trCxB) + " = 0", p4Y, 2, trBLUE, "dSSE/db = 0");
    TEXTOUT_P("Solved exactly (CAS solve), beside the decimals:", G1, 10, p4Y, 1, #000000h, 302);
    p4Y := p4Y + 11;
    p4Y := trCxFrac("b", trCxShow(trCxSb), trB, p4Y);
    p4Y := trCxFrac("a", trCxShow(trCxSa), trA0, p4Y);
  ELSE
    // ---- decimal: why, then the normal equations with the sums ----
    p4S := "The CAS did not answer: no exact fractions.";
    IF trCxSt == 3 THEN p4S := "A value needs over 6 decimals: no exact fractions."; END;
    IF trCxSt == 4 THEN p4S := "Too long to send to the CAS here."; END;
    TEXTOUT_P(p4S, G1, 10, p4Y, 1, trRED, 302);
    p4Y := p4Y + 11;
    p4Sx := 0; p4Sxx := 0; p4Sxa := 0;
    FOR p4I FROM 1 TO trK DO
      p4Sx := p4Sx + trTot(p4I) * trX(p4I);
      p4Sxx := p4Sxx + trTot(p4I) * trX(p4I)^2;
      p4Sxa := p4Sxa + trX(p4I) * trA(p4I);
    END;
    TEXTOUT_P("Both slopes of the SSE = 0, the normal equations:", G1, 10, p4Y, 1, #000000h, 302);
    TEXTOUT_P("N a + (sum n x) b = A", G1, 10, p4Y + 11, 1, trBLUE, 302);
    TEXTOUT_P("  " + trNum(trNT) + " a + " + trNum(p4Sx) + " b = " + trNum(trAT), G1, 10, p4Y + 22, 1, trBLUE, 302);
    TEXTOUT_P("(sum n x) a + (sum n x^2) b = sum x a", G1, 10, p4Y + 33, 1, trBLUE, 302);
    TEXTOUT_P("  " + trNum(p4Sx) + " a + " + trNum(p4Sxx) + " b = " + trNum(p4Sxa), G1, 10, p4Y + 44, 1, trBLUE, 302);
    TEXTOUT_P("Solving them (the program's decimals):", G1, 10, p4Y + 55, 1, #000000h, 302);
    TEXTOUT_P("  b = " + trDec(trB) + ",  a = " + trDec(trA0), G1, 10, p4Y + 66, 1, trMAG, 302);
    p4Y := p4Y + 77;
  END;
  // ---- what the line means (at most y 189: the footer is at 205) ----
  TEXTOUT_P("= least squares on all N people's 0/1 outcomes;", G1, 10, p4Y, 1, trGRAY, 302);
  TEXTOUT_P("the trend " + trChi + " asks whether this b is 0.", G1, 10, p4Y + 11, 1, trGRAY, 302);
  IF trCxSt == 1 AND DIM(trCxNt) > 0 AND p4Y + 22 <= 189 THEN
    TEXTOUT_P("Exact: score " + trCxNt + "; nothing rounded.", G1, 10, p4Y + 22, 1, trGRAY, 302);
  END;
  TEXTOUT_P("4/4, ⇆ or tap = page, Esc = menu", G1, 10, 205, 1, trFOOT, 300);
END;


// ------------------------------------------------------------
// PLOT 1: proportions by level, Wilson bars, WLS trend line and
// its 95% band; a dot can be dragged
// ------------------------------------------------------------

// pixel x of a score, for the plot window xlo..xhi
trPxOf(pxV, pxLo, pxHi)
BEGIN
  RETURN trPL + IP((pxV - pxLo) / (pxHi - pxLo) * (trPR - trPL) + 0.5);
END;

// pixel y of a proportion on plot 1's scale (0 .. trYhi), kept
// inside the plot box -- during a drag the scale is frozen, so a
// Wilson bar or the band can run past the top and is cut there
trPyOf(pyV)
BEGIN
  LOCAL pyY;
  pyY := trPB - IP(pyV / trYhi * (trPB - trPTop) + 0.5);
  RETURN MAX(trPTop, MIN(trPB, pyY));
END;

// plot 1's y-scale: a nice step, the top a whole number of steps
// above the largest Wilson upper limit (never above 1)
trFitY()
BEGIN
  trYStp := trNice(trHiMx / 4);
  trYhi := CEILING(trHiMx / trYStp - 0.000001) * trYStp;
  IF trYhi > 1 THEN trYhi := 1; END;
  trYDec := MAX(0, -FLOOR(LOG(trYStp)));
END;

// Compose plot 1 in G1.  dpF is the fraction of the fitted slope to
// draw: 1 normally; 0 .. 1 while the fit animates (trAnimate).
trDrawProp(dpF)
BEGIN
  LOCAL dpI, dpX, dpY, dpY1, dpY2, dpV;
  LOCAL dpCol, dpS, dpPx, dpXv, dpN, dpB;
  LOCAL dpA, dpFit, dpH, dpSe, dpR, dpT;
  LOCAL dpV1, dpV2, dpTa, dpTb, dpLo, dpHi;
  // ---- the y window: refitted unless a drag has frozen it ----
  IF trFrz == 0 THEN trFitY(); END;
  // the line drawn in this frame: slope trB * dpF, pivoting on the
  // point (xbar, overall p) that every weighted line passes through
  dpB := trB * dpF;
  dpA := trAT / trNT - dpB * trXbar;

  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  IF trIllu == 1 THEN
    TEXTOUT_P("Made-up data: proportion by level", G1, 10, 4, 2, trRED);
  ELSE
    TEXTOUT_P("Proportion with the outcome, by level", G1, 10, 4, 2, #000000h);
  END;
  // ---- grid and y ticks ----
  dpN := 0;
  dpV := 0;
  WHILE dpV <= trYhi + 0.000001 AND dpN < 12 DO
    dpY := trPyOf(dpV);
    LINE_P(G1, trPL, dpY, trPR, dpY, trGRID);
    TEXTOUT_P(STRING(dpV, 2, trYDec), G1, 6, dpY - 5, 1, trGRAY, 38);
    dpV := dpV + trYStp;
    dpN := dpN + 1;
  END;
  // ---- the 95% band round the fitted line, light magenta, one
  // pixel column at a time (only for the fitted line, not while it
  // animates): fit +/- 1.96 * sqrt(P(1-P) (1/N + (x-xbar)^2/Sxx)).
  // Narrowest at xbar, where the line is pinned by all N people. ----
  IF dpF >= 1 THEN
    FOR dpPx FROM trPL + 1 TO trPR DO
      dpXv := trXlo + (dpPx - trPL) / (trPR - trPL) * (trXhi - trXlo);
      dpFit := trA0 + trB * dpXv;
      dpH := trZ * (trPQ * (1 / trNT + (dpXv - trXbar)^2 / trSxx))^0.5;
      dpY1 := trPyOf(dpFit + dpH);
      dpY2 := trPyOf(dpFit - dpH);
      IF dpY2 > dpY1 THEN LINE_P(G1, dpPx, dpY1, dpPx, dpY2, trLMAG); END;
    END;
  END;
  // ---- axes ----
  LINE_P(G1, trPL, trPTop, trPL, trPB, #000000h);
  LINE_P(G1, trPL, trPB, trPR, trPB, #000000h);
  // ---- x ticks: one per level, labelled with its score ----
  FOR dpI FROM 1 TO trK DO
    dpX := trPxOf(trX(dpI), trXlo, trXhi);
    LINE_P(G1, dpX, trPB, dpX, trPB + 3, #000000h);
    TEXTOUT_P(trXStr(trX(dpI)), G1, dpX - 8, trPB + 5, 1, trGRAY, 40);
  END;
  TEXTOUT_P("score:", G1, 6, trPB + 5, 1, trGRAY);
  // ---- overall proportion, dotted ----
  dpY := trPyOf(trAT / trNT);
  trDotLine(G1, trPL, trPR, dpY, trGRAY);
  TEXTOUT_P("overall", G1, trPR - 44, dpY - 11, 1, trGRAY);
  // ---- the line, 2 px thick: one LINE_P between the points where
  // it enters and leaves the window 0..trYhi.  Along the box, t runs
  // 0 (left edge, value dpV1) to 1 (right edge, value dpV2). ----
  dpV1 := dpA + dpB * trXlo;
  dpV2 := dpA + dpB * trXhi;
  dpLo := 0; dpHi := 1;
  IF dpV2 == dpV1 THEN
    IF dpV1 < 0 OR dpV1 > trYhi THEN dpHi := -1; END;    // flat, outside
  ELSE
    dpTa := (0 - dpV1) / (dpV2 - dpV1);                  // where it crosses 0
    dpTb := (trYhi - dpV1) / (dpV2 - dpV1);              // ... and the top
    dpLo := MAX(0, MIN(dpTa, dpTb));
    dpHi := MIN(1, MAX(dpTa, dpTb));
  END;
  IF dpHi > dpLo THEN
    dpX := trPL + IP(dpLo * (trPR - trPL) + 0.5);
    dpT := trPL + IP(dpHi * (trPR - trPL) + 0.5);
    dpY1 := trPyOf(dpV1 + dpLo * (dpV2 - dpV1));
    dpY2 := trPyOf(dpV1 + dpHi * (dpV2 - dpV1));
    LINE_P(G1, dpX, dpY1, dpT, dpY2, trMAG);
    LINE_P(G1, dpX, MAX(trPTop, dpY1 - 1), dpT, MAX(trPTop, dpY2 - 1), trMAG);
  END;
  // ---- Wilson bars and the dots, coloured by the residual from
  // this frame's line: blue above it by more than 1 SE, red below
  // by more than 1 SE, grey within; SE = sqrt(P(1-P)/n).  The dots'
  // pixel centres go into trDX, trDY for the touch handler. ----
  trDX := {}; trDY := {};
  FOR dpI FROM 1 TO trK DO
    dpX := trPxOf(trX(dpI), trXlo, trXhi);
    dpY1 := trPyOf(trHi(dpI));
    dpY2 := trPyOf(trLo(dpI));
    dpY := trPyOf(trP(dpI));
    dpR := trP(dpI) - (dpA + dpB * trX(dpI));
    dpSe := (trPQ / trTot(dpI))^0.5;
    dpCol := trGRAY;
    IF dpR > dpSe THEN dpCol := trBLUE; END;
    IF dpR < -dpSe THEN dpCol := trRED; END;
    LINE_P(G1, dpX, dpY1, dpX, dpY2, dpCol);
    LINE_P(G1, dpX - 3, dpY1, dpX + 3, dpY1, dpCol);
    LINE_P(G1, dpX - 3, dpY2, dpX + 3, dpY2, dpCol);
    trDot(G1, dpX, dpY, dpCol);
    IF dpI == trSel THEN                         // the selected level: a black ring
      ARC_P(G1, dpX, dpY, 6, #000000h);
      ARC_P(G1, dpX, dpY, 7, #000000h);
    END;
    trDX := CONCAT(trDX, {dpX});
    trDY := CONCAT(trDY, {dpY});
  END;
  // ---- legend ----
  TEXTOUT_P("above", G1, trPL + 4, trPTop + 2, 1, trBLUE);
  TEXTOUT_P("below", G1, trPL + 38, trPTop + 2, 1, trRED);
  TEXTOUT_P("within 1 SE", G1, trPL + 72, trPTop + 2, 1, trGRAY);
  TEXTOUT_P("WLS line, 95% band", G1, trPL + 146, trPTop + 2, 1, trMAG);
  // ---- footer, four lines above the soft keys (y 170-217) ----
  // the selected level: its count, p, the line's value there, and
  // the residual in SE units
  dpFit := dpA + dpB * trX(trSel);
  dpR := (trP(trSel) - dpFit) / (trPQ / trTot(trSel))^0.5;
  IF dpR < 0 THEN dpS := "-"; ELSE dpS := "+"; END;
  // beta 3c: trFx (ROUND, then fixed decimals), so a small value never
  // prints in scientific form -- the emulator showed "(+9.6E-2 SE)"
  dpS := "L" + STRING(trSel, 2, 0) + ": x = " + trXStr(trX(trSel)) + ", " + STRING(trA(trSel), 2, 0) + "/" + STRING(trTot(trSel), 2, 0) + " = " + trFx(trP(trSel), 3) + ", line " + trFx(dpFit, 3) + " (" + dpS + trFx(ABS(dpR), 1) + " SE)";
  TEXTOUT_P(dpS, G1, 6, 170, 1, #000000h, 310);
  // the trend statistic: while the fit animates, the value it would
  // have for the slope drawn (it grows with the slope squared)
  IF dpF >= 1 THEN
    TEXTOUT_P(trChi + "trend = " + trFx(trX2T, 2) + ", " + trPStr(trPT) + ", scores " + trScName(), G1, 6, 182, 1, trMAG, 310);
  ELSE
    TEXTOUT_P("fitting: slope " + CHAR(215) + " " + STRING(dpF, 2, 2) + ", " + trChi + "trend = " + trFx(trX2T * dpF^2, 2), G1, 6, 182, 1, trMAG, 310);
  END;
  TEXTOUT_P("Pearson " + trChi + " = " + trFx(trX2P, 2), G1, 6, 194, 1, trBLUE, 120);
  IF trDFD >= 1 THEN
    TEXTOUT_P("departure = " + trFx(trX2D, 2) + ", " + trPStr(trPD), G1, 130, 194, 1, trGRN, 186);
  ELSE
    TEXTOUT_P("departure: needs 3+ levels", G1, 130, 194, 1, trGRN, 186);
  END;
  TEXTOUT_P("drag a dot or ⇅ = cases, ⇆ = level, Esc = menu", G1, 6, 206, 1, trFOOT, 310);
END;

// The fit, animated: the line starts flat at the overall proportion
// and swings to the fitted slope over 15 frames while the trend
// statistic counts up.  Any key or tap skips to the end.
trAnimate()
BEGIN
  LOCAL anI, anGo, anK, anM, anM1;
  anI := 0; anGo := 1;
  WHILE anGo == 1 DO
    trDrawProp(anI / 15);
    trShow();
    anK := GETKEY;
    anM := MOUSE;
    anM1 := anM(1);
    IF anK > -1 OR SIZE(anM1) > 0 THEN anGo := 0; END;
    anI := anI + 1;
    IF anI >= 15 THEN anGo := 0; END;
  END;
  REPEAT anK := GETKEY; UNTIL anK == -1;         // a skipping key or
  trDrain();                                     // tap does nothing else
END;


// ------------------------------------------------------------
// PLOT 2: odds ratios vs level 1, log scale, Woolf bars
// ------------------------------------------------------------
trDrawOR()
BEGIN
  LOCAL doI, doLlo, doLhi;
  LOCAL doSpan, doX, doY, doY1, doY2, doCol;
  LOCAL doD, doM, doTk, doS, doV, doMul;
  // ---- window: x as in plot 1; y = log10(OR) ----
  doLlo := 0; doLhi := 0;                        // always include OR = 1
  FOR doI FROM 2 TO trK DO
    doLlo := MIN(doLlo, LOG(trORL(doI)));
    doLhi := MAX(doLhi, LOG(trORH(doI)));
  END;
  doSpan := MAX(0.5, doLhi - doLlo);
  doLlo := doLlo - 0.08 * doSpan; doLhi := doLhi + 0.08 * doSpan;

  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  IF trIllu == 1 THEN
    TEXTOUT_P("Made-up data: OR vs level 1 (log)", G1, 10, 4, 2, trRED);
  ELSE
    TEXTOUT_P("Odds ratio vs level 1 (log scale)", G1, 10, 4, 2, #000000h);
  END;
  // ---- ticks at 1, 2, 5 in each decade ----
  doMul := {1, 2, 5};
  FOR doD FROM FLOOR(doLlo) TO CEILING(doLhi) DO
    FOR doM FROM 1 TO 3 DO
      doTk := doD + LOG(doMul(doM));
      IF doTk >= doLlo AND doTk <= doLhi THEN
        doY := trPB - IP((doTk - doLlo) / (doLhi - doLlo) * (trPB - trPTop) + 0.5);
        LINE_P(G1, trPL, doY, trPR, doY, trGRID);
        doV := doMul(doM) * 10^doD;
        IF doV >= 1 THEN
          TEXTOUT_P(STRING(doV, 2, 0), G1, 6, doY - 5, 1, trGRAY, 38);
        ELSE
          TEXTOUT_P(STRING(doV, 2, MAX(1, -doD)), G1, 6, doY - 5, 1, trGRAY, 38);
        END;
      END;
    END;
  END;
  // ---- axes and the OR = 1 line ----
  LINE_P(G1, trPL, trPTop, trPL, trPB, #000000h);
  LINE_P(G1, trPL, trPB, trPR, trPB, #000000h);
  doY := trPB - IP((0 - doLlo) / (doLhi - doLlo) * (trPB - trPTop) + 0.5);
  trDotLine(G1, trPL, trPR, doY, trGRAY);
  TEXTOUT_P("OR = 1", G1, trPR - 40, doY - 11, 1, trGRAY);
  // ---- x ticks ----
  FOR doI FROM 1 TO trK DO
    doX := trPxOf(trX(doI), trXlo, trXhi);
    LINE_P(G1, doX, trPB, doX, trPB + 3, #000000h);
    TEXTOUT_P(trXStr(trX(doI)), G1, doX - 8, trPB + 5, 1, trGRAY, 40);
  END;
  TEXTOUT_P("score:", G1, 6, trPB + 5, 1, trGRAY);
  // ---- level 1: a hollow circle at OR = 1 ----
  doX := trPxOf(trX(1), trXlo, trXhi);
  doCol := trGRAY;
  IF trSel == 1 THEN doCol := trRED; END;
  ARC_P(G1, doX, doY, 3, doCol);
  TEXTOUT_P("ref", G1, doX - 8, doY + 6, 1, doCol);
  // ---- the other levels: Woolf bars and dots ----
  FOR doI FROM 2 TO trK DO
    doX := trPxOf(trX(doI), trXlo, trXhi);
    doY := trPB - IP((LOG(trOR(doI)) - doLlo) / (doLhi - doLlo) * (trPB - trPTop) + 0.5);
    doY1 := trPB - IP((LOG(trORH(doI)) - doLlo) / (doLhi - doLlo) * (trPB - trPTop) + 0.5);
    doY2 := trPB - IP((LOG(trORL(doI)) - doLlo) / (doLhi - doLlo) * (trPB - trPTop) + 0.5);
    doCol := trMAG;
    IF doI == trSel THEN doCol := trRED; END;
    LINE_P(G1, doX, doY1, doX, doY2, doCol);
    LINE_P(G1, doX - 3, doY1, doX + 3, doY1, doCol);
    LINE_P(G1, doX - 3, doY2, doX + 3, doY2, doCol);
    trDot(G1, doX, doY, doCol);
  END;
  // ---- footer: readout, note, why the axis is logarithmic, keys ----
  IF trSel == 1 THEN
    doS := "L1: the reference level, OR = 1";
  ELSE
    doS := "L" + STRING(trSel, 2, 0) + " vs L1:  OR " + trOrStr(trOR(trSel)) + "  (" + trOrStr(trORL(trSel)) + " - " + trOrStr(trORH(trSel)) + ")";
  END;
  TEXTOUT_P(doS, G1, 6, 170, 1, trRED, 310);
  IF trHald == 1 THEN
    TEXTOUT_P("+0.5 in every cell (zero cell); Woolf 95% CI", G1, 6, 182, 1, trRED, 310);
  ELSE
    TEXTOUT_P("Woolf 95% CI", G1, 6, 182, 1, trGRAY, 310);
  END;
  TEXTOUT_P("Log scale: OR 0.5 and 2 are equally far from 1.", G1, 6, 194, 1, #000000h, 310);
  TEXTOUT_P("⇆ or tap = level, Esc = menu", G1, 6, 206, 1, trFOOT, 310);
END;


// ------------------------------------------------------------
// TOUCH ON THE PLOTS
// ------------------------------------------------------------

// A finger came down at (ptX0, ptY0) on plot ptW (4 proportions,
// 5 odds ratios).  On plot 1, within 10 px of a dot: that dot is
// grabbed -- if the finger then moves more than 8 px it is a drag,
// and every sample sets the level's cases from the finger's height
// on the FROZEN y-scale; if it lifts without moving it was a tap,
// which selects the level.  Anywhere else in the plot box, a tap
// selects the nearest level in x (within 14 px).
trPlotTouch(ptW, ptX0, ptY0)
BEGIN
  LOCAL ptI, ptHit, ptD, ptX, ptM, ptM1;
  LOCAL ptOn, ptMv, ptNew, ptCLo, ptCHi, ptP;
  ptHit := 0;
  IF ptW == 4 THEN
    ptD := 11;
    FOR ptI FROM 1 TO trK DO
      ptX := MAX(ABS(ptX0 - trDX(ptI)), ABS(ptY0 - trDY(ptI)));
      IF ptX <= 10 AND ptX < ptD THEN ptD := ptX; ptHit := ptI; END;
    END;
  END;
  IF ptHit > 0 THEN
    trSel := ptHit;
    // the other levels' cases stay fixed; keep at least one case and
    // one non-case in the whole table, or V = 0 and the test has
    // nothing to compare
    ptCLo := MAX(0, 1 - (trAT - trA(ptHit)));
    ptCHi := MIN(trTot(ptHit), trNT - 1 - (trAT - trA(ptHit)));
    ptMv := 0;
    trFrz := 1;                                  // freeze the y-scale
    REPEAT
      ptM := MOUSE;
      ptM1 := ptM(1);
      ptOn := 0;                                 // finger still on the screen?
      IF SIZE(ptM1) > 0 THEN
        IF ptM1(1) >= 0 AND ptM1(1) <= 319 AND ptM1(2) >= 0 AND ptM1(2) <= 239 THEN ptOn := 1; END;
      END;
      IF ptOn == 1 THEN
        IF ABS(ptM1(1) - ptX0) > 8 OR ABS(ptM1(2) - ptY0) > 8 THEN ptMv := 1; END;
        IF ptMv == 1 THEN
          // finger height -> proportion on the frozen scale -> cases
          ptP := (trPB - ptM1(2)) / (trPB - trPTop) * trYhi;
          ptP := MAX(0, MIN(trYhi, ptP));
          ptNew := MAX(ptCLo, MIN(ptCHi, ROUND(ptP * trTot(ptHit), 0)));
          IF ptNew <> trA(ptHit) THEN            // redraw only on a change
            trA(ptHit) := ptNew;
            trEdit := 1;
            trCompute();
            trDrawProp(1);
            trShow();
          END;
        END;
      END;
    UNTIL ptOn == 0;
    trFrz := 0;                                  // refit on lift
  ELSE
    ptD := 15;
    IF ptY0 >= trPTop AND ptY0 <= trPB THEN
      FOR ptI FROM 1 TO trK DO
        ptX := trPxOf(trX(ptI), trXlo, trXhi);
        IF ABS(ptX0 - ptX) < ptD THEN ptD := ABS(ptX0 - ptX); ptHit := ptI; END;
      END;
    END;
    IF ptHit > 0 THEN trSel := ptHit; END;
  END;
  trDrain();
END;


// ------------------------------------------------------------
// THE VIEW LOOP: results pages 1-4 and the two plots, one loop,
// one soft-key strip.  vwV: 1-3 results pages, 4 proportions plot,
// 5 odds-ratio plot, 6 results page 4 (least squares as calculus,
// beta 3 -- numbered 6 so that the plots keep 4 and 5 and the drag
// code is untouched).
//   Soft keys: Prop, OR, Scores, Data, Restore, Menu.
//   Keys:  results  Left/Right = page, Enter = next (menu after 4)
//          plots    Left/Right = level, Up/Down = cases (plot 1),
//                   Enter = menu, - = the other plot
//          all      + = Scores, 0 = Restore, Esc = menu
//   Taps:  results  left half = back, right half = on
//          plots    a dot = grab and drag it (plot 1); else select
// ------------------------------------------------------------

// The results pages in the order the student turns them: 1, 2, 3,
// then 6 (page 4).  trPgStep(v, d) is the page d steps (+1 or -1)
// from v, staying at the first and the last.
trPgStep(psV, psD)
BEGIN
  LOCAL psO, psJ, psK;
  psO := {1, 2, 3, 6};
  psK := 1;
  FOR psJ FROM 1 TO 4 DO
    IF psO(psJ) == psV THEN psK := psJ; END;
  END;
  psK := MAX(1, MIN(4, psK + psD));
  RETURN psO(psK);
END;

trDrawView(dwV)
BEGIN
  CASE
    IF dwV == 1 THEN trPage1(); END;
    IF dwV == 2 THEN trPageChi(); END;
    IF dwV == 3 THEN trPage3(); END;
    IF dwV == 4 THEN trDrawProp(1); END;
    IF dwV == 6 THEN
      // page 4: the CAS is asked here only if the table or the scores
      // changed since it last answered (trCxMake keeps its answers)
      trCxMake();
      trPage4();
    END;
    DEFAULT trDrawOR();
  END;
  trShow();
END;

trView(vwStart)
BEGIN
  LOCAL vwV, vwOld, vwGo, vwK, vwM, vwM1;
  LOCAL vwSl, vwStp, vwLo, vwHi, vwA;
  vwV := vwStart; vwOld := 0; vwGo := 1;
  trFrz := 0;
  REPEAT vwK := GETKEY; UNTIL vwK == -1;
  trDrain();
  WHILE vwGo == 1 DO
    // entering the proportions plot from anywhere else: animate the fit
    IF vwV == 4 AND vwOld <> 4 THEN trAnimate(); END;
    vwOld := vwV;
    trDrawView(vwV);
    // wait for a key or a touch (poll both)
    REPEAT
      vwK := GETKEY;
      vwM := MOUSE;
      vwM1 := vwM(1);
    UNTIL vwK > -1 OR SIZE(vwM1) > 0;
    vwSl := 0;                                   // the soft key to act on
    IF SIZE(vwM1) > 0 THEN
      IF vwM1(2) >= 218 THEN
        vwSl := MIN(6, IP(vwM1(1) / 53.4) + 1);  // the soft-key strip
        trDrain();
      ELSE
        IF vwV <= 3 OR vwV == 6 THEN
          // a results page: left half = back, right half = on
          IF vwM1(1) < 160 THEN vwV := trPgStep(vwV, -1); ELSE vwV := trPgStep(vwV, 1); END;
          trDrain();
        ELSE
          trPlotTouch(vwV, vwM1(1), vwM1(2));    // drains on lift
        END;
      END;
    ELSE
      IF vwK == 4 THEN vwGo := 0; END;           // Esc = menu
      IF vwK == 50 THEN vwSl := 3; END;          // +   = Scores
      IF vwK == 42 THEN vwSl := 5; END;          // 0   = Restore
      IF vwV <= 3 OR vwV == 6 THEN
        IF vwK == 7 THEN vwV := trPgStep(vwV, -1); END;
        IF vwK == 8 THEN vwV := trPgStep(vwV, 1); END;
        IF vwK == 30 THEN
          IF vwV == 6 THEN vwGo := 0; ELSE vwV := trPgStep(vwV, 1); END;
        END;
        IF vwK == 45 THEN vwV := 4; END;         // -   = the plots
      ELSE
        IF vwK == 7 THEN trSel := MAX(1, trSel - 1); END;
        IF vwK == 8 THEN trSel := MIN(trK, trSel + 1); END;
        IF vwK == 30 THEN vwGo := 0; END;
        // Up/Down on plot 1: the keyboard's drag.  One press moves the
        // dot about 1/40 of the y-scale (at least one case).
        IF vwV == 4 THEN
          IF vwK == 2 OR vwK == 12 THEN
            vwStp := MAX(1, ROUND(trTot(trSel) * trYhi / 40, 0));
            vwLo := MAX(0, 1 - (trAT - trA(trSel)));
            vwHi := MIN(trTot(trSel), trNT - 1 - (trAT - trA(trSel)));
            vwA := trA(trSel);
            IF vwK == 2 THEN vwA := MIN(vwHi, vwA + vwStp); END;
            IF vwK == 12 THEN vwA := MAX(vwLo, vwA - vwStp); END;
            IF vwA <> trA(trSel) THEN
              trA(trSel) := vwA;
              trEdit := 1;
              trCompute();
            END;
          END;
        END;
        IF vwK == 45 THEN vwV := 9 - vwV; END;   // -   = the other plot
      END;
    END;
    // ---- the soft keys (and their hardware shortcuts) ----
    IF vwSl == 1 THEN vwV := 4; END;             // Prop
    IF vwSl == 2 THEN vwV := 5; END;             // OR
    IF vwSl == 3 THEN                            // Scores: next set, recompute
      trScM := trScM + 1;
      IF trScM > 2 THEN trScM := 0; END;
      trSetX();
      trCompute();
    END;
    IF vwSl == 4 THEN                            // Data: the entry dialogs
      IF trEnter() == 1 THEN vwV := 1; END;
      // the tap that closed a dialog can arrive again on the screen
      // underneath, so drain keys and touch before the next frame
      REPEAT vwK := GETKEY; UNTIL vwK == -1;
      trDrain();
    END;
    IF vwSl == 5 THEN                            // Restore the data as entered
      trA := trCopy(trOA);
      trScM := 0;
      trSetX();
      trEdit := 0;
      trCompute();
    END;
    IF vwSl == 6 THEN vwGo := 0; END;            // Menu
  END;
  REPEAT vwK := GETKEY; UNTIL vwK == -1;
END;


// ------------------------------------------------------------
// HELP: three pages
// ------------------------------------------------------------
trHelp()
BEGIN
  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 235, trBLUE, trLBLU);
  TEXTOUT_P("Chi-square for trend: the idea", G1, 12, 8, 2, trBLUE);
  TEXTOUT_P("Exposure levels have an ORDER (dose, hours,", G1, 12, 32, 1, #000000h);
  TEXTOUT_P("grade).  Pearson's chi-square (k-1 df) asks", G1, 12, 45, 1, #000000h);
  TEXTOUT_P("only 'do the proportions differ?'.  The trend", G1, 12, 58, 1, #000000h);
  TEXTOUT_P("test (1 df) asks 'do they rise or fall", G1, 12, 71, 1, #000000h);
  TEXTOUT_P("steadily with the score?'", G1, 12, 84, 1, #000000h);
  TEXTOUT_P("1 df concentrates the evidence: more power", G1, 12, 103, 1, trMAG);
  TEXTOUT_P("than Pearson when the line really fits.", G1, 12, 116, 1, trMAG);
  TEXTOUT_P("Pearson - trend = departure from linearity", G1, 12, 129, 1, trGRN);
  TEXTOUT_P("(k-2 df); large means 'not a straight line'.", G1, 12, 142, 1, trGRN);
  TEXTOUT_P("Scores: use the doses themselves when known;", G1, 12, 161, 1, trGRAY);
  TEXTOUT_P("otherwise 1, 2, 3 ...  The trend answer depends", G1, 12, 174, 1, trGRAY);
  TEXTOUT_P("on the spacing; Pearson's does not.  Press", G1, 12, 187, 1, trGRAY);
  TEXTOUT_P("Scores: 1,2,3 / as entered / midranks.", G1, 12, 200, 1, trGRAY);
  TEXTOUT_P("any key = page 2, Esc = menu", G1, 12, 220, 1, trFOOT);
  BLIT_P(G0, 0, 0, G1);
  // beta 3c: Esc leaves help from any page (the emulator turned the page)
  IF trWaitAny() == 4 THEN RETURN; END;

  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 235, trBLUE, trLBLU);
  TEXTOUT_P("Reading the plots", G1, 12, 8, 2, trBLUE);
  TEXTOUT_P("Proportions plot: a dot per level, p = cases/n", G1, 12, 28, 1, trBLUE);
  TEXTOUT_P("with its 95% Wilson bar; the weighted (WLS)", G1, 12, 41, 1, trBLUE);
  TEXTOUT_P("line in magenta with its 95% band; overall p", G1, 12, 54, 1, trBLUE);
  TEXTOUT_P("dotted.  Dots: blue above the line, red below,", G1, 12, 67, 1, trBLUE);
  TEXTOUT_P("grey within 1 SE of it (SE^2 = P(1-P)/n).", G1, 12, 80, 1, trBLUE);
  TEXTOUT_P("Drag a dot up or down (or Up/Down) to change", G1, 12, 96, 1, #000000h);
  TEXTOUT_P("its cases: the line tilts, the band and the", G1, 12, 109, 1, #000000h);
  TEXTOUT_P("tests follow.  Restore puts the data back.", G1, 12, 122, 1, #000000h);
  TEXTOUT_P("Odds ratio plot: each level against level 1,", G1, 12, 138, 1, trMAG);
  TEXTOUT_P("log scale: 0.5 and 2 lie equally far from 1;", G1, 12, 151, 1, trMAG);
  TEXTOUT_P("Woolf 95% intervals; dotted line = OR 1.", G1, 12, 164, 1, trMAG);
  TEXTOUT_P("Band: fit +/- 1.96 SE, SE^2 = P(1-P) times", G1, 12, 180, 1, trGRAY);
  TEXTOUT_P("(1/N + (x - xbar)^2 / sum n(x - xbar)^2).", G1, 12, 193, 1, trGRAY);
  TEXTOUT_P("Results page 4: where the line comes from.", G1, 12, 206, 1, trMAG);
  TEXTOUT_P("any key = page 3, Esc = menu", G1, 12, 220, 1, trFOOT);
  BLIT_P(G0, 0, 0, G1);
  IF trWaitAny() == 4 THEN RETURN; END;

  DIMGROB_P(G1, 320, 240, #FFFFFFh);
  RECT_P(G1, 4, 4, 315, 235, trBLUE, trLBLU);
  TEXTOUT_P("Keys, and the same test elsewhere", G1, 12, 8, 2, trBLUE);
  TEXTOUT_P("Soft keys: Prop and OR = the two plots;", G1, 12, 28, 1, #000000h);
  TEXTOUT_P("Scores = 1,2,3 / as entered / midranks; Data", G1, 12, 41, 1, #000000h);
  TEXTOUT_P("= type new counts; Restore = undo the drags;", G1, 12, 54, 1, #000000h);
  TEXTOUT_P("Menu = back.  Keys: ⇆ = level, ⇅ = cases,", G1, 12, 67, 1, #000000h);
  TEXTOUT_P("+ = Scores, 0 = Restore, - = other plot,", G1, 12, 80, 1, #000000h);
  TEXTOUT_P("Esc = menu.  Results: 4 pages, ⇆ or tap.", G1, 12, 93, 1, #000000h);
  TEXTOUT_P("The proportions plot animates the fit as it", G1, 12, 106, 1, trMAG);
  TEXTOUT_P("opens (any key or tap skips).  Page 4: the CAS", G1, 12, 119, 1, trMAG);
  TEXTOUT_P("finds the line's slope as an exact fraction.", G1, 12, 132, 1, trMAG);
  // beta 3d: these five lines now say WHICH form each program prints
  // (fleet review 19 Oct 2026, C-T2): Epi Info prints Mantel's
  // (N-1) r^2, as this program does; R and SAS print Armitage's N r^2,
  // 72.688 for the snoring table against 72.659 here.  Mantel's form
  // stays the headline.  JMP is no longer named (its form was not
  // checked).
  TEXTOUT_P("Same test elsewhere: Epi Info StatCalc prints", G1, 12, 146, 1, trGRAY);
  TEXTOUT_P("Mantel's (1963) (N-1) r^2, as this page does;", G1, 12, 159, 1, trGRAY);
  TEXTOUT_P("R prop.trend.test and SAS PROC FREQ / TREND", G1, 12, 172, 1, trGRAY);
  TEXTOUT_P("print Armitage's (1955) N r^2: 72.688 for", G1, 12, 185, 1, trGRAY);
  TEXTOUT_P("snoring, 72.659 here.  The ratio is N/(N-1).", G1, 12, 198, 1, trGRAY);
  TEXTOUT_P("any key = menu", G1, 12, 220, 1, trFOOT);
  BLIT_P(G0, 0, 0, G1);
  trWaitAny();
END;


// ------------------------------------------------------------
// INIT AND MAIN
// ------------------------------------------------------------
trInit()
BEGIN
  trBLUE := #003CC8h; trMAG := #C800C8h; trRED := #FF0000h;
  trGRN := #008C3Ch; trORNG := #FF4500h; trGRAY := #5A5A5Ah;
  trFOOT := #1E90FFh; trGRID := #E1E4EBh; trLORN := #FFF0DCh;
  trLBLU := #EBF2FFh; trLMAG := #FAEBFAh;
  // chi (U+03C7) + superscript two (U+00B2), as in ChiTable
  trChi := CHAR(967) + CHAR(178);
  // plot box: the footer's four lines (170-217) and the soft keys'
  // strip (218-239) sit below it
  trPL := 46; trPR := 306; trPTop := 28; trPB := 152;
  trZ := 1.959964;
  trHas := 0; trSel := 1; trK := 0;
  trX := {}; trA := {}; trTot := {}; trLbl := ""; trLbl2 := "";
  trOX := {}; trOA := {}; trScM := 0; trEdit := 0; trIllu := 0;
  trYhi := 1; trYStp := 0.2; trYDec := 1; trFrz := 0;
  trDX := {}; trDY := {};
  // beta 3: least squares as calculus
  trCasOK := 0; trCasMx := 4000; trCxK := ""; trCxSt := 0; trCxNt := "";
  trCxE := ""; trCxA := ""; trCxB := ""; trCxSa := ""; trCxSb := "";
END;

EXPORT TrendTest()
BEGIN
  LOCAL mnI, mnRun, mnS;
  trInit();
  mnRun := 1;
  WHILE mnRun == 1 DO
    mnI := 0;
    // beta 3d: "Plot: proportions + trend line" (30 characters) was at
    // the edge of the list-form box, which cut a 35-character line
    // after about 31 on the emulator; now 29 (fleet review 19 Oct
    // 2026, C-T3).  Comments go above the CHOOSE, never inside the
    // list (a Syntax Error on the Prime).
    CHOOSE(mnI, "Chi-square for linear trend", {
      "Enter data (2-6 levels)",
      "Results",
      "Plot: proportions, trend line",
      "Plot: odds ratios vs level 1",
      "Examples",
      "Help",
      "Quit"});
    CASE
      IF mnI == 0 OR mnI == 7 THEN mnRun := 0; END;
      IF mnI == 1 THEN
        IF trEnter() == 1 THEN trView(1); END;
      END;
      IF mnI >= 2 AND mnI <= 4 AND trHas == 0 THEN
        MSGBOX("No data yet.  Use Enter data, or load one of the Examples, first.");
      END;
      IF mnI == 2 AND trHas == 1 THEN trView(1); END;
      IF mnI == 3 AND trHas == 1 THEN trView(4); END;
      IF mnI == 4 AND trHas == 1 THEN trView(5); END;
      IF mnI == 5 THEN
        IF trExamples() == 1 THEN trView(1); END;
      END;
      IF mnI == 6 THEN trHelp(); END;
    END;
  END;
  // Hand the summary back as the program's value: Home shows it in
  // the history and no Terminal screen opens.  The numbers stay in
  // TR_X2T TR_P TR_X2 TR_DF TR_OR TR_B.
  IF trHas == 1 THEN
    mnS := "TrendTest: X2 trend = " + trFx(trX2T, 4) + " (1 df), " + trPStr(trPT) + "; Pearson X2 = " + trFx(trX2P, 4) + " (" + STRING(trK - 1, 2, 0) + " df)";
    IF trScM <> 0 THEN mnS := mnS + "; scores " + trScName(); END;
    IF trEdit == 1 THEN mnS := mnS + "; counts edited"; END;
    RETURN mnS;
  END;
  RETURN "TrendTest: no table";
END;
