# imports section
import glob
import os
import re
import warnings

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
from tcrdist.adpt_funcs import import_adaptive_file

warnings.filterwarnings("ignore")

SEED = 42                 # controls rarefaction draws, permutations and the bootstrap
N_RAREFY = 20             # rarefaction replicates averaged per sample
MIN_TB_CLONOTYPES = 5     # below this, TB-specific clonality is not estimable
N_PERM = 10000            # permutations for the stratified test
N_BOOT = 2000             # bootstrap resamples for the effect-size CI

SAMPLE_DIR = "datasets/experimental/samples"
CACHE_DIR = "cache/samples_imgt"
GROUP_ORDER = ['Progressor', 'Controller']

# auxiliary functions necessary for analysis
def standardize_gene_name(series):
    '''Converts V and J genes names to proper IMGT format'''
    series = series.astype(str).str.strip().str.upper().str.split('*').str[0]
    series = series.str.replace(r'TRBV0', 'TRBV', regex = True)
    series = series.str.replace(r'TRBJ0', 'TRBJ', regex = True)
    series = series.str.replace(r'TCRBV', 'TRBV', regex = True)
    series = series.str.replace(r'TCRBJ', 'TRBJ', regex = True)
    series = series.str.replace(r'(TRB[V])(\d+)$', r'\1\2-1', regex=True)
    series = series.str.replace(r'(TRB[J])(\d+)$', r'\1\2-1', regex=True)
    return series

def standardize_CDR3(series):
    '''Standardizes the CDR3s as well for good measure'''
    return series.astype(str).str.strip().str.upper()

def canonical_sample_id(s):
    '''Collapses the two sample-ID spellings onto one join key.

    The sample key writes ACS IDs as `04-0333-D0`; the exported files are named
    `04-0333_D0`. Treating '_' and '-' as the same separator makes the join
    total instead of silently discarding most of the ACS cohort.'''
    return re.sub(r'[_\-]+', '-', str(s).strip().upper())

def shannon_clonality(counts):
    '''Normalized Shannon clonality, 1 - H/ln(N). NaN below two clonotypes.

    With a single clonotype the normalisation is 0/0. The first version returned
    1.0 there, which invents a maximally clonal sample out of one we could not
    actually measure, and those fabricated 1.0s sat at the top of the figure.'''
    c = np.asarray(counts, dtype=float)
    c = c[c > 0]
    n = c.size
    if n < 2:
        return np.nan
    p = c / c.sum()
    entropy = -np.sum(p * np.log(p))
    # A perfectly even subset gives entropy == ln(N) exactly; floating point can
    # push that a hair past 1.0, so clamp rather than emit a tiny negative.
    return max(0.0, 1.0 - entropy / np.log(n))

def rarefied_sample_metrics(counts, is_tb, depth, n_rep, rng):
    '''Depth-standardized repertoire statistics for one sample.

    Both clonality and the frequency of database-matched clonotypes increase
    with sequencing depth, because a deeper library simply observes more of the
    repertoire. Drawing a fixed number of templates without replacement -- a
    multivariate hypergeometric draw from the clonotype urn -- puts every sample
    on the same footing; averaging over replicates removes draw-to-draw noise.'''
    counts = np.asarray(counts, dtype=np.int64)
    is_tb = np.asarray(is_tb, dtype=bool)
    if counts.sum() < depth:
        return np.nan, np.nan, np.nan, np.nan
    tot_clon, tb_freq, tb_clon, tb_n = [], [], [], []
    for _ in range(n_rep):
        drawn = rng.multivariate_hypergeometric(counts, depth, method='marginals')
        tb_drawn = drawn[is_tb]
        tot_clon.append(shannon_clonality(drawn))
        tb_freq.append(tb_drawn.sum() / depth)
        tb_n.append(int((tb_drawn > 0).sum()))
        if (tb_drawn > 0).sum() >= MIN_TB_CLONOTYPES:
            tb_clon.append(shannon_clonality(tb_drawn))

    # A sample only gets a TB-specific clonality if most of its draws actually
    # cleared the threshold. Averaging over just the replicates that happened to
    # reach it selects the luckiest draws and hands back exactly the unstable,
    # inflated values this metric is prone to -- the same pathology as scoring a
    # one-clonotype subset as perfectly clonal.
    tb_clonality = (float(np.mean(tb_clon))
                    if len(tb_clon) >= n_rep / 2 else np.nan)
    return (float(np.mean(tot_clon)), float(np.mean(tb_freq)),
            tb_clonality, float(np.mean(tb_n)))

def format_p_value(p):
    '''Formatting the P value string from the permutation test'''
    if np.isnan(p):
        return 'P n/a'
    if p < 0.001:
        return 'P < 0.001'
    return f'P = {p:.3f}'

# stratified (matched-set) inference
def stratified_rank_stat(values, is_case, strata):
    '''Van Elteren style stratified Wilcoxon statistic.

    Values are ranked *within* each matched set and the standardized rank-sums
    are added up, so a progressor is only ever compared against the controllers
    they were matched to and between-set variation never enters the contrast.'''
    total = 0.0
    for s in np.unique(strata):
        m = strata == s
        v, c = values[m], is_case[m]
        n, n1 = v.size, int(c.sum())
        if n1 == 0 or n1 == n or n < 2:
            continue
        r = stats.rankdata(v)
        exp = n1 * (n + 1) / 2.0
        var = n1 * (n - n1) * (n + 1) / 12.0
        if var <= 0:
            continue
        total += (r[c].sum() - exp) / np.sqrt(var)
    return total

def stratified_perm_test(values, is_case, strata, n_perm=N_PERM, seed=SEED):
    '''Two-sided p-value, permuting group labels within each matched set.

    Permuting inside the strata is what makes this a test of the matched design:
    the null is "progressor status is exchangeable within a matched set", and no
    distributional assumption is needed.

    Within a stratum the ranks depend only on the values, so they are computed
    once and every permutation just re-selects which of them the cases hold.
    That lets all n_perm permutations for a stratum be drawn as one array.'''
    ok = ~np.isnan(values)
    values, is_case, strata = values[ok], is_case[ok], strata[ok]
    informative = np.array([is_case[strata == s].sum() > 0 and (~is_case[strata == s]).sum() > 0
                            for s in strata])
    values, is_case, strata = values[informative], is_case[informative], strata[informative]
    if values.size == 0:
        return np.nan, np.nan, 0, 0

    rng = np.random.default_rng(seed)
    obs = 0.0
    perm_totals = np.zeros(n_perm)
    for s in np.unique(strata):
        m = strata == s
        v, c = values[m], is_case[m]
        n, n1 = v.size, int(c.sum())
        if n1 == 0 or n1 == n or n < 2:
            continue
        r = stats.rankdata(v)
        exp = n1 * (n + 1) / 2.0
        sd = np.sqrt(n1 * (n - n1) * (n + 1) / 12.0)
        if sd <= 0:
            continue
        obs += (r[c].sum() - exp) / sd
        # One random arrangement of this stratum's labels per permutation.
        picks = np.argsort(rng.random((n_perm, n)), axis = 1)[:, :n1]
        perm_totals += (r[picks].sum(axis = 1) - exp) / sd

    ge = int((np.abs(perm_totals) >= abs(obs)).sum())
    # (r + 1) / (n + 1): a finite permutation set can never justify p = 0.
    return obs, (ge + 1) / (n_perm + 1), int(is_case.sum()), int((~is_case).sum())

def hodges_lehmann(values, is_case, strata, n_boot=N_BOOT, seed=SEED):
    '''Within-set Hodges-Lehmann shift with a bootstrap CI over matched sets.

    The median of every progressor-minus-controller difference formed inside a
    matched set. A negative result is only informative with an interval attached:
    the CI is what says how large a difference the data can still rule out.'''
    ok = ~np.isnan(values)
    values, is_case, strata = values[ok], is_case[ok], strata[ok]

    # Each matched set's progressor-minus-controller differences, computed once.
    per_set = {}
    for s in np.unique(strata):
        m = strata == s
        a, b = values[m & is_case], values[m & ~is_case]
        if a.size and b.size:
            per_set[s] = (a[:, None] - b[None, :]).ravel()
    if not per_set:
        return np.nan, np.nan, np.nan

    sets = np.array(list(per_set))
    point = float(np.median(np.concatenate([per_set[s] for s in sets])))

    # Resampling whole matched sets keeps the bootstrap at the level the design
    # randomises at, rather than treating individual donors as exchangeable.
    rng = np.random.default_rng(seed)
    boots = np.empty(n_boot)
    for i in range(n_boot):
        pick = rng.choice(sets, size = sets.size, replace = True)
        boots[i] = np.median(np.concatenate([per_set[s] for s in pick]))
    lo, hi = np.percentile(boots, [2.5, 97.5])
    return point, float(lo), float(hi)

def bh_fdr(pvals):
    '''Benjamini-Hochberg adjusted p-values across the three metrics.'''
    p = np.asarray(pvals, dtype=float)
    ok = ~np.isnan(p)
    q = np.full(p.shape, np.nan)
    pp = p[ok]
    n = pp.size
    order = np.argsort(pp)
    ranked = pp[order] * n / (np.arange(n) + 1)
    ranked = np.minimum.accumulate(ranked[::-1])[::-1]
    out = np.empty(n)
    out[order] = np.clip(ranked, 0, 1)
    q[ok] = out
    return q

# external reference set collection/cleaning/standardization/combining
iedb = pd.read_csv("datasets/reference/iedb.csv")
vdjdb = pd.read_csv("datasets/reference/vdjdb.tsv", sep = '\t')

iedb.dropna(inplace = True)
iedb = iedb[iedb['Epitope - Source Organism'].str.contains("Mycobacterium tuberculosis")]
iedb["V"] = standardize_gene_name(iedb["Chain 2 - Curated V Gene"])
iedb["J"] = standardize_gene_name(iedb["Chain 2 - Curated J Gene"])
iedb["CDR3"] = standardize_CDR3(iedb["Chain 2 - CDR3 Curated"])
iedb = iedb[['CDR3', 'V', 'J']]

vdjdb = vdjdb[['CDR3', 'V', 'J']]
vdjdb.dropna(inplace = True)
vdjdb["V"] = standardize_gene_name(vdjdb["V"])
vdjdb["J"] = standardize_gene_name(vdjdb["J"])
vdjdb["CDR3"] = standardize_CDR3(vdjdb["CDR3"])

ext_set = pd.concat([vdjdb, iedb], ignore_index=True)
ext_set = ext_set[ext_set['V'].str.startswith('TRB', na = False)].copy()
ext_set = ext_set[ext_set['J'].str.startswith('TRB', na = False)].copy()
ext_set.drop_duplicates(keep = 'first', inplace = True)

ref_keys = set(map(tuple, ext_set[['V', 'J', 'CDR3']].to_numpy()))
print(f"reference clonotypes (IEDB Mtb + VDJdb, exact V/J/CDR3): {len(ref_keys):,}")

reference clonotypes (IEDB Mtb + VDJdb, exact V/J/CDR3): 1,426

# one-time parse of the raw Adaptive exports into a compact cache
#
# The 277 raw .tsv exports total ~110 GB. import_adaptive_file (the same loader
# the first version used) is applied once per file and the result -- collapsed
# onto the (V, J, CDR3-aa) clonotype -- is cached to ~120 MB, so re-running the
# analysis costs seconds instead of re-reading 110 GB.
os.makedirs(CACHE_DIR, exist_ok=True)

def cache_one(path):
    name = os.path.basename(path)[:-9]              # strip '_TCRB.tsv'
    out = os.path.join(CACHE_DIR, f"{name}.tsv.gz")
    if os.path.exists(out):
        return out
    df = import_adaptive_file(path)
    df = df.dropna(subset = ['cdr3_b_aa', 'v_b_gene', 'j_b_gene', 'templates'])
    df = df[~df['cdr3_b_aa'].str.contains(r'\*', na = False)]
    # Collapse nucleotide-level rearrangements onto the amino-acid clonotype.
    # The reference match is defined at (V, J, CDR3-aa), so that is the unit the
    # analysis should count; leaving them separate inflates the clonotype count
    # and therefore deflates normalised clonality.
    g = (df.groupby(['v_b_gene', 'j_b_gene', 'cdr3_b_aa'], observed = True)
           .agg(templates = ('templates', 'sum'), n_rearr = ('templates', 'size'))
           .reset_index())
    g.to_csv(out, sep = '\t', index = False, compression = 'gzip')
    return out

raw_files = sorted(glob.glob(os.path.join(SAMPLE_DIR, '*.tsv')))
for i, f in enumerate(raw_files, 1):
    cache_one(f)
    if i % 50 == 0:
        print(f"  cached {i}/{len(raw_files)}")
print(f"cache ready: {len(glob.glob(os.path.join(CACHE_DIR, '*.tsv.gz')))} samples")

  cached 50/277
  cached 100/277
  cached 150/277
  cached 200/277
  cached 250/277
cache ready: 277 samples

# join the sample key to the files, reporting every ID that fails to match
meta = pd.read_csv("datasets/experimental/metadata/samplekey.csv")
cached = {canonical_sample_id(os.path.basename(f)[:-7]): f
          for f in glob.glob(os.path.join(CACHE_DIR, '*.tsv.gz'))}

meta['canon'] = meta['Sample.ID'].map(canonical_sample_id)
meta['path'] = meta['canon'].map(cached)

key_without_file = meta.loc[meta['path'].isna(), 'Sample.ID'].tolist()
file_without_key = sorted(set(cached) - set(meta['canon']))

print(f"sample-key records: {len(meta)}; joined to a file: {meta['path'].notna().sum()}")
print(f"key rows with no file ({len(key_without_file)}): {key_without_file}")
print(f"files with no key row ({len(file_without_key)}): {file_without_key}")

meta = meta.dropna(subset = ['path']).copy()

sample-key records: 272; joined to a file: 271
key rows with no file (1): ['8_BL']
files with no key row (6): ['1552-M18', '2864-M18', '2864-M6', '3094-M6', '8-BL-1', '907-BL-TC']

# primary analysis set: one baseline sample per donor, inside matched sets
#
# Fix 1 (pseudoreplication): 67 donors contributed 2-3 longitudinal samples to
# the first version's pooled test. Restricting to the baseline visit makes every
# row an independent person.
# Fix 2 (the matched design): each progressor was matched to controllers on
# time-to-diagnosis within a cohort, so (Cohort, Days.To.TB) identifies the
# matched set that the stratified test conditions on.
base = meta[meta['VISIT'].isin(['Month0', 'Day0'])].copy()
base = base.sort_values('Sample.ID').drop_duplicates('Donor.ID', keep = 'first')
base['set_id'] = base['Cohort'].astype(str) + '|' + base['Days.To.TB'].astype(str)

n_sets = base['set_id'].nunique()
usable = (base.groupby('set_id')['Group']
              .agg(lambda g: g.eq('Progressor').any() and g.eq('Controller').any()))
print(f"baseline donors: {len(base)}  "
      f"({(base.Group == 'Progressor').sum()} progressor, "
      f"{(base.Group == 'Controller').sum()} controller)")
print(f"donors contributing more than one row: {int(base['Donor.ID'].duplicated().sum())}")
print(f"matched sets: {n_sets}; sets with both arms present: {int(usable.sum())}")

# Rarefaction depth: the 5th percentile of library size, so nearly every sample
# is retained while all of them are compared at one common sequencing depth.
depths = [int(pd.read_csv(p, sep = '\t', usecols = ['templates'])['templates'].sum())
          for p in base['path']]
base['total_templates'] = depths
DEPTH = int(np.percentile(depths, 5))
print(f"\nlibrary size: min={min(depths):,}  median={int(np.median(depths)):,}  "
      f"max={max(depths):,}  (a {max(depths)/min(depths):.1f}x spread)")
print(f"rarefying every sample to {DEPTH:,} templates")

baseline donors: 140  (52 progressor, 88 controller)
donors contributing more than one row: 0
matched sets: 51; sets with both arms present: 42

library size: min=2,925  median=238,463  max=530,352  (a 181.3x spread)
rarefying every sample to 43,414 templates

# per-sample metrics, all computed on the rarefied repertoire
rng = np.random.default_rng(SEED)
rows = []

for _, r in base.iterrows():
    d = pd.read_csv(r['path'], sep = '\t')
    V = standardize_gene_name(d['v_b_gene'])
    J = standardize_gene_name(d['j_b_gene'])
    C = standardize_CDR3(d['cdr3_b_aa'])
    is_tb = np.fromiter((k in ref_keys for k in zip(V, J, C)), bool, len(d))

    tot, freq, tb_clon, tb_n = rarefied_sample_metrics(
        d['templates'].to_numpy(), is_tb, DEPTH, N_RAREFY, rng)

    rows.append({'Sample.ID': r['Sample.ID'], 'Donor.ID': r['Donor.ID'],
                 'Group': r['Group'], 'set_id': r['set_id'],
                 'Age': r['Age'], 'Sex': r['Sex'], 'Cohort': r['Cohort'],
                 'total_templates': r['total_templates'], 'n_clonotypes': len(d),
                 'tb_clonotypes_obs': int(is_tb.sum()),
                 'Total Clonality': tot, 'Total TB Frequency': freq,
                 'TB Clonality': tb_clon, 'tb_clonotypes_rarefied': tb_n})

res = pd.DataFrame(rows)
# Keep the original figure's group order rather than order of appearance.
res['Group'] = pd.Categorical(res['Group'], categories = GROUP_ORDER, ordered = True)

os.makedirs('results', exist_ok = True)
res.to_csv('results/per_sample_metrics.tsv', sep = '\t', index = False)

print(f"samples below the rarefaction depth (excluded): "
      f"{int(res['Total Clonality'].isna().sum())}")
print(f"TB-matched clonotypes per sample (before rarefaction): "
      f"median={res['tb_clonotypes_obs'].median():.0f}, "
      f"IQR {res['tb_clonotypes_obs'].quantile(.25):.0f}-"
      f"{res['tb_clonotypes_obs'].quantile(.75):.0f}, "
      f"max={res['tb_clonotypes_obs'].max():.0f}")
print(f"TB clonality estimable (>= {MIN_TB_CLONOTYPES} matched clonotypes in the "
      f"majority of draws): {int(res['TB Clonality'].notna().sum())} of {len(res)} samples")
res.head()

samples below the rarefaction depth (excluded): 7
TB-matched clonotypes per sample (before rarefaction): median=16, IQR 10-26, max=50
TB clonality estimable (>= 5 matched clonotypes in the majority of draws): 93 of 140 samples

# stratified testing, effect sizes and multiplicity adjustment
METRICS = ['Total TB Frequency', 'TB Clonality', 'Total Clonality']

case = (res['Group'] == 'Progressor').to_numpy()
strata = res['set_id'].to_numpy()
summary = []

for col in METRICS:
    v = res[col].to_numpy(float)
    stat, p_strat, n1, n0 = stratified_perm_test(v, case, strata)
    hl, lo, hi = hodges_lehmann(v, case, strata)
    ok = ~np.isnan(v)
    p_mw = stats.mannwhitneyu(v[ok & case], v[ok & ~case],
                              alternative = 'two-sided').pvalue
    summary.append({'Metric': col, 'n prog': n1, 'n cont': n0,
                    'median prog': np.nanmedian(v[case]),
                    'median cont': np.nanmedian(v[~case]),
                    'HL shift': hl, 'CI low': lo, 'CI high': hi,
                    'P (stratified)': p_strat, 'P (pooled MW)': p_mw})

summary = pd.DataFrame(summary)
summary['Q (BH)'] = bh_fdr(summary['P (stratified)'].to_numpy())
summary.to_csv('results/summary_statistics.tsv', sep = '\t', index = False)

pd.set_option('display.width', 200)
print(summary.to_string(index = False, float_format = lambda x: f'{x:.4g}'))

            Metric  n prog  n cont  median prog  median cont  HL shift     CI low   CI high  P (stratified)  P (pooled MW)  Q (BH)
Total TB Frequency      44      73    0.0001647    0.0001866 8.062e-06 -3.183e-05 5.183e-05          0.3961         0.8901  0.6607
      TB Clonality      28      45      0.07781      0.04715  0.001853    -0.1445    0.1453          0.6607         0.9396  0.6607
   Total Clonality      44      73       0.0759       0.0764  0.001993    -0.0228   0.02145          0.5236         0.3646  0.6607

# compiles data and plots figure
plot_df = res.melt(id_vars = ['Group'], value_vars = METRICS,
                   var_name = 'Metric', value_name = 'Value').dropna(subset = ['Value'])

fig, axes = plt.subplots(1, 3, figsize = (12, 6), sharey = False)
palette = {'Progressor': 'firebrick', 'Controller': 'steelblue'}
titles = {'Total TB Frequency': 'Total TB-Specific Clone Frequency',
          'TB Clonality': 'Clonality of TB-Specific Subset',
          'Total Clonality': 'Clonality of Total Sample'}
ylabels = {'Total TB Frequency': 'Total Frequency (rarefied)',
           'TB Clonality': 'Normalized Clonality (0: Diverse, 1: Monoclonal)',
           'Total Clonality': 'Normalized Clonality (0: Diverse, 1: Monoclonal)'}

for ax, col in zip(axes, METRICS):
    sub = plot_df[plot_df['Metric'] == col]
    sns.boxplot(data = sub, x = 'Group', y = 'Value', hue = 'Group', legend = False,
                order = GROUP_ORDER, hue_order = GROUP_ORDER,
                palette = palette, width = 0.5, ax = ax,
                medianprops = {'color': 'black'})
    row = summary[summary['Metric'] == col].iloc[0]
    ax.set_title(titles[col])
    ax.set_ylabel(ylabels[col])
    # Headroom for the annotation, so it never overlaps a whisker or outlier.
    lo, hi = ax.get_ylim()
    ax.set_ylim(lo, hi + 0.42 * (hi - lo))
    ax.text(0.5, 0.88,
            f"Stratified Wilcoxon:\n{format_p_value(row['P (stratified)'])}"
            f"  (q = {row['Q (BH)']:.3f})\n"
            f"n = {int(row['n prog'])} vs {int(row['n cont'])}",
            transform = ax.transAxes, ha = 'center', fontsize = 10,
            bbox = dict(boxstyle = 'round', facecolor = 'white', alpha = 0.85))

fig.suptitle('Comparison of TB-Specific/Nonspecific TCR Repertoires: '
             'Progressors vs. Controllers\n'
             'baseline sample per donor, rarefied to a common depth, '
             'tested within matched sets', fontsize = 13)
plt.tight_layout()
plt.show()

#	Issue in the first version	What this version does
1	Pseudoreplication. 67 donors contributed 2-3 longitudinal samples, all treated as independent observations.	One baseline sample per donor, so every observation is an independent person.
2	The matched design was discarded. Each progressor was matched to controllers on time-to-diagnosis, age band and sex; the pooled test ignores that.	A stratified (van Elteren) rank test, permuting group labels within each matched set.
3	Clonality is confounded with sequencing depth. Normalised Shannon clonality divides by `ln(N)`, and deeper libraries observe more rare clonotypes. Library size spans a 181-fold range across these donors.	Every sample is rarefied to a common template count before any statistic is computed.
4	A group-assignment bug. In the TB-specific loop, the progressor branch appended its missing value to the controller list.	Removed; missing values are reported, never reassigned.
5	Missing values were imputed to the group median, which shrinks variance toward the null.	Missing values are excluded and counted, not filled.
6	A one-clonotype subset returned clonality = 1.0, inventing a maximally clonal sample from an unmeasurable one.	Returns `NaN` below two clonotypes; TB-specific clonality additionally requires ≥5 matched clonotypes in the majority of rarefaction draws, so the estimate cannot rest on a few lucky ones.
7	71 of 277 sample files silently failed to join to the sample key, dropping most of the ACS cohort. The two sources spell IDs as `04-0333_D0` and `04-0333-D0`.	Separator-insensitive join, and every unmatched ID on either side is printed.
8	No effect sizes, no CIs, no multiplicity adjustment across the three tests.	Hodges-Lehmann shift with a bootstrap CI, plus Benjamini-Hochberg across the three metrics.

	Sample.ID	Donor.ID	Group	set_id	Age	Sex	Cohort	total_templates	n_clonotypes	tb_clonotypes_obs	Total Clonality	Total TB Frequency	TB Clonality	tb_clonotypes_rarefied
0	01-0935-D0	01-0935	Controller	ACS\|633	13	female	ACS	49412	33596	11	0.080004	0.000359	0.058767	10.00
1	04-0333-D0	04-0333	Controller	ACS\|256	17	female	ACS	101945	65188	12	0.120644	0.000242	0.146504	6.20
2	04-0695-D0	04-0695	Progressor	ACS\|617	14	male	ACS	77451	32710	2	0.249073	0.000038	NaN	1.35
3	04-0699-D0	04-0699	Controller	ACS\|557	14	female	ACS	75664	34580	5	0.211949	0.055852	NaN	3.70
4	04-0741-D0	04-0741	Controller	ACS\|211	14	female	ACS	138524	79214	13	0.148387	0.000236	0.067163	6.30

Metric	n (prog vs cont)	Median prog	Median cont	Hodges-Lehmann shift (95% CI)	P (stratified)	Q (BH)	P (pooled MW)
Total TB Frequency	44 vs 73	0.0001647	0.0001866	+8.06e-06 (-3.18e-05 to +5.18e-05)	0.396	0.661	0.890
TB Clonality	28 vs 45	0.07781	0.04715	+0.00185 (-0.144 to +0.145)	0.661	0.661	0.940
Total Clonality	44 vs 73	0.0759	0.0764	+0.00199 (-0.0228 to +0.0215)	0.524	0.661	0.365

TCR Analysis: TB Progression/Control Clonality¶

Naraen Palanikumar¶

Revision note — what changed and why¶

Conclusions¶