feat: 导入 IRSTD 研究代码与文档（小目标检测）

- 新增 irstd 核心包、scripts、tests、configs
- 新增研究文档：literature_review / mathematical_analysis / novelty_matrix / experiment_plan / research_decision / research_log / manuscript_outline
- 新增 evidence、results 数据与审计产物
- 新增 pyproject.toml、requirements、.gitignore

irstd/__init__.py

@@ -0,0 +1 @@
+"""Research prototype; no pretrained weights or benchmark results included."""

irstd/core.py

@@ -0,0 +1,229 @@
+"""Reference CPU implementation. Coordinates are (row, col), native pixels."""
+import itertools
+import math
+import numpy as np
+from scipy import ndimage
+from scipy.stats import beta, binom
+
+
+def components(mask):
+    mask = np.asarray(mask, dtype=bool)
+    if mask.ndim != 2:
+        raise ValueError('Expected a 2D mask')
+    labels, n = ndimage.label(mask, structure=np.ones((3, 3), dtype=int))
+    return [labels == i for i in range(1, n + 1)]
+
+
+def centers(masks):
+    return np.asarray([np.argwhere(m).mean(0) for m in masks], float).reshape(-1, 2)
+
+
+def eligibility(pred_masks, gt_masks, radius=3.0):
+    if radius < 0:
+        raise ValueError('radius must be nonnegative')
+    a, b = centers(pred_masks), centers(gt_masks)
+    return np.linalg.norm(a[:, None] - b[None, :], axis=-1) <= radius
+
+
+def matching(edges):
+    """Deterministic maximum cardinality bipartite matching, not greedy matching."""
+    edges = np.asarray(edges, bool)
+    if edges.ndim != 2:
+        raise ValueError('Expected two-dimensional eligibility')
+    owner = [-1] * edges.shape[1]
+    def augment(i, seen):
+        for j in np.flatnonzero(edges[i]):
+            if j in seen:
+                continue
+            seen.add(j)
+            if owner[j] < 0 or augment(owner[j], seen):
+                owner[j] = i
+                return True
+        return False
+    for i in range(edges.shape[0]):
+        augment(i, set())
+    return sorted((i, j) for j, i in enumerate(owner) if i >= 0)
+
+
+def candidate_pool(prob, base_threshold=0.1, cap=64):
+    """Threshold ONCE; subsequent gate thresholds delete whole fixed components."""
+    prob = np.asarray(prob, float)
+    if prob.ndim != 2 or not np.isfinite(prob).all() or np.any((prob < 0) | (prob > 1)):
+        raise ValueError('Expected finite 2D probabilities in [0,1]')
+    if not 0 <= base_threshold <= 1 or not isinstance(cap, int) or cap < 1:
+        raise ValueError('Invalid threshold/cap')
+    masks = components(prob >= base_threshold)
+    masks.sort(key=lambda m: -float(prob[m].max()))  # stable ties: raster component order
+    masks = masks[:cap]
+    features = []
+    for m in masks:
+        p = prob[m]
+        ring = ndimage.binary_dilation(m, iterations=2) & ~m
+        contrast = float(p.mean() - prob[ring].mean()) if ring.any() else 0.0
+        features.append([float(p.max()), float(p.mean()), float(p.std()),
+                         float(np.log1p(m.sum())), contrast])
+    return masks, np.asarray(features, np.float32).reshape(-1, 5)
+
+
+def canonical_labels(pred_masks, gt_mask, radius=3.0):
+    gt = components(gt_mask)
+    pairs = matching(eligibility(pred_masks, gt, radius))
+    labels = np.zeros(len(pred_masks), np.float32)
+    areas = np.ones(len(pred_masks), np.float32)
+    for i, j in pairs:
+        labels[i] = 1
+        areas[i] = gt[j].sum()
+    return labels, areas
+
+
+def image_stats(pred_masks, gt_mask, radius=3.0):
+    gt_mask = np.asarray(gt_mask, bool)
+    gt = components(gt_mask)
+    pairs = matching(eligibility(pred_masks, gt, radius))
+    tp = len(pairs)
+    union_pred = np.zeros_like(gt_mask)
+    for m in pred_masks:
+        if m.shape != gt_mask.shape:
+            raise ValueError('Shape mismatch')
+        union_pred |= m
+    intersection = int((union_pred & gt_mask).sum())
+    union = int((union_pred | gt_mask).sum())
+    fp_pixels = int((union_pred & ~gt_mask).sum())
+    matched = {i for i, _ in pairs}
+    # Separate legacy-style unmatched-component area from literal background-pixel FP.
+    unmatched_area = sum(int(m.sum()) for i, m in enumerate(pred_masks) if i not in matched)
+    return dict(tp=tp, fp=len(pred_masks)-tp, fn=len(gt)-tp,
+                intersection=intersection, union=union,
+                pred_pixels=int(union_pred.sum()), gt_pixels=int(gt_mask.sum()),
+                fp_pixels=fp_pixels, unmatched_area=unmatched_area,
+                pixels=gt_mask.size, niou=(intersection/union if union else 1.0))
+
+
+def aggregate(rows):
+    if not rows:
+        raise ValueError('Cannot evaluate an empty split')
+    sums = {k: sum(r[k] for r in rows) for k in rows[0] if k != 'niou'}
+    def ratio(a, b):
+        return a/b if b else None
+    tp, fp, fn = (sums[k] for k in ('tp', 'fp', 'fn'))
+    return dict(images=len(rows), **sums,
+                IoU=ratio(sums['intersection'], sums['union']),
+                nIoU=float(np.mean([r['niou'] for r in rows])),
+                Pd=ratio(tp, tp+fn), object_F1=ratio(2*tp, 2*tp+fp+fn),
+                pixel_F1=ratio(2*sums['intersection'], sums['pred_pixels']+sums['gt_pixels']),
+                FP_per_image=fp/len(rows), P_any_FA=sum(r['fp']>0 for r in rows)/len(rows),
+                Fa_pixel=sums['fp_pixels']/sums['pixels'],
+                Fa_unmatched_area=sums['unmatched_area']/sums['pixels'])
+
+
+def cp_upper(successes, n, error):
+    if n < 1 or not 0 < error < 1 or not 0 <= successes <= n:
+        raise ValueError('Invalid binomial parameters')
+    return 1.0 if successes == n else float(beta.ppf(1-error, successes+1, n-successes))
+
+
+def calibrate(fp_counts, thresholds, cap, alpha=0.1, delta=0.05, risk='count', method='bonferroni'):
+    """Finite-grid CP calibration. Rows MUST be iid held-out deployment units.
+
+    Does not infer iid sampling from filenames. Empty fallback is explicit.
+    """
+    counts = np.asarray(fp_counts)
+    ts = np.asarray(thresholds, float)
+    if (counts.ndim != 2 or counts.shape[0] < 1 or counts.shape[1] != len(ts)
+        or len(ts) < 1 or not np.isfinite(counts).all()
+        or not np.isfinite(ts).all() or np.any((ts<0)|(ts>1))
+        or len(np.unique(ts)) != len(ts) or not isinstance(cap, int) or cap < 1
+        or np.any(counts < 0) or np.any(counts > cap) or np.any(counts != np.floor(counts))
+        or not 0 < delta < 1 or not 0 <= alpha <= (cap if risk=='count' else 1)
+        or risk not in ('count', 'any') or method not in ('bonferroni', 'fixed_sequence')):
+        raise ValueError('Invalid calibration input')
+    n, m = counts.shape
+    if method == 'fixed_sequence':
+        # Order fixed before inspecting labels; stop at the first failed test.
+        order = np.argsort(-ts).tolist()
+        tested = []
+        chosen = None
+        limit = cap if risk == 'count' else 1
+        for j in order:
+            values = counts[:, j] if risk == 'count' else (counts[:, j] > 0)
+            mean = float(np.mean(values))/limit
+            target = alpha/limit
+            if target == 1:
+                pvalue = 0.0  # Risk <= known support is deterministic.
+            elif target == 0 or mean >= target:
+                pvalue = 1.0
+            elif risk == 'any':
+                pvalue = float(binom.cdf(int(np.sum(values)), n, target))
+            else:
+                kl = ((mean*math.log(mean/target) if mean else 0.0)
+                      + (1-mean)*math.log((1-mean)/(1-target)))
+                pvalue = math.exp(-n*kl)
+            tested.append(dict(threshold=float(ts[j]), pvalue=pvalue))
+            if pvalue > delta:
+                break
+            chosen = j
+        return dict(method=method, risk=risk, alpha=alpha, delta=delta, n=n, cap=cap,
+                    thresholds=ts.tolist(), tested=tested,
+                    threshold=None if chosen is None else float(ts[chosen]),
+                    reject_all=chosen is None,
+                    certified_risk_budget=0.0 if chosen is None else alpha,
+                    status='reject_all_no_nontrivial_certificate' if chosen is None else 'certified_under_iid_assumptions')
+    tails = cap if risk == 'count' else 1
+    bounds = [sum(cp_upper(int((counts[:,j]>=k).sum()), n, delta/(m*tails))
+                  for k in range(1, tails+1)) for j in range(m)]
+    feasible = [j for j in range(m) if bounds[j] <= alpha]
+    # Lowest threshold keeps largest nested set; no label-based utility search needed.
+    chosen = min(feasible, key=lambda j: ts[j]) if feasible else None
+    return dict(method=method, risk=risk, alpha=alpha, delta=delta, n=n, cap=cap,
+                thresholds=ts.tolist(), upper_bounds=bounds,
+                threshold=None if chosen is None else float(ts[chosen]),
+                reject_all=chosen is None,
+                selected_upper_bound=0.0 if chosen is None else bounds[chosen],
+                status='reject_all_no_nontrivial_certificate' if chosen is None else 'certified_under_iid_assumptions')
+
+
+def softplus_bound(scores, negative, tau=1.0):
+    if tau <= 0:
+        raise ValueError('tau must be positive')
+    return float(np.logaddexp(0, np.asarray(scores)[np.asarray(negative, bool)]/tau).sum()/math.log(2))
+
+
+def exact_mwis(weights, edges):
+    """Tiny-instance oracle only; exponential, capped at 20 nodes."""
+    w = np.asarray(weights, float)
+    if len(w) > 20:
+        raise ValueError('Oracle limited to 20 nodes')
+    best, selection = 0.0, []
+    for bits in itertools.product((0,1), repeat=len(w)):
+        if any(bits[i] and bits[j] for i,j in edges):
+            continue
+        score = float(w @ bits)
+        if score > best:
+            best, selection = score, [i for i,b in enumerate(bits) if b]
+    return best, selection
+
+
+def best_path(node_scores, transition_costs):
+    """Complete path through layered DAG. inf cost denotes forbidden edge."""
+    if not node_scores or len(transition_costs) != len(node_scores)-1:
+        raise ValueError('Inconsistent layers')
+    dp = np.asarray(node_scores[0], float)
+    if not len(dp):
+        return -math.inf, []
+    parents = []
+    for scores, costs in zip(node_scores[1:], transition_costs):
+        scores, costs = np.asarray(scores,float), np.asarray(costs,float)
+        if costs.shape != (len(dp), len(scores)) or not len(scores):
+            raise ValueError('Invalid transitions')
+        val = dp[:,None]-costs
+        par = val.argmax(0)
+        dp = scores + val[par, np.arange(len(scores))]
+        parents.append(par)
+    k = int(dp.argmax())
+    best = float(dp[k])
+    if not np.isfinite(best):
+        return -math.inf, []
+    path = [k]
+    for par in reversed(parents):
+        k = int(par[k]); path.append(k)
+    return best, path[::-1]

irstd/data.py

@@ -0,0 +1,44 @@
+import csv
+import hashlib
+from pathlib import Path
+import numpy as np
+from PIL import Image
+
+def sha(path):
+    return hashlib.sha256(Path(path).read_bytes()).hexdigest()
+
+def read_manifest(path):
+    path=Path(path)
+    with path.open(newline='') as stream:
+        rows=list(csv.DictReader(stream))
+    if not rows or not {'id','image','mask','split','group'} <= set(rows[0]):
+        raise ValueError('Manifest requires id,image,mask,split,group')
+    ids=set(); groups={}; hashes={}
+    for r in rows:
+        if r['id'] in ids or not r['id'] or '/' in r['id'] or '\\' in r['id'] or r['id'] in ('.','..'):
+            raise ValueError('Duplicate/unsafe id')
+        ids.add(r['id'])
+        if r['split'] not in ('train','dev','cal','test') or not r['group']:
+            raise ValueError('Invalid split/group')
+        if r['group'] in groups and groups[r['group']]!=r['split']:
+            raise ValueError('Scene group crosses splits')
+        groups[r['group']]=r['split']
+        for k in ('image','mask'):
+            r[k]=str((path.parent/r[k]).resolve())
+            if not Path(r[k]).is_file(): raise FileNotFoundError(r[k])
+        h=sha(r['image'])
+        if h in hashes and hashes[h]!=r['split']:
+            raise ValueError('Identical image bytes cross splits')
+        hashes[h]=r['split']
+    return rows
+
+def read_pair(row):
+    with Image.open(row['image']) as img:
+        image=np.asarray(img).copy()
+    with Image.open(row['mask']) as img:
+        mask=np.asarray(img).copy()
+    if image.ndim!=2 or image.dtype!=np.uint8:
+        raise ValueError('Prototype requires native-resolution 8-bit grayscale; no implicit conversion')
+    if mask.ndim!=2 or mask.shape!=image.shape:
+        raise ValueError('Mask shape/format mismatch')
+    return image.astype(np.float32)/255., mask>0

irstd/models.py

@@ -0,0 +1,51 @@
+"""Original compact baseline, NOT a reproduction of MSHNet or DNANet."""
+import torch
+from torch import nn
+from torch.nn import functional as F
+
+class TinyUNet(nn.Module):
+    def __init__(self, width=16):
+        super().__init__()
+        def block(a,b):
+            return nn.Sequential(nn.Conv2d(a,b,3,padding=1),nn.ReLU(),
+                                 nn.Conv2d(b,b,3,padding=1),nn.ReLU())
+        self.enc = block(1,width)
+        self.mid = block(width,width*2)
+        self.dec = block(width*3,width)
+        self.head = nn.Conv2d(width,1,1)
+    def forward(self,x):
+        a = self.enc(x)
+        b = self.mid(F.max_pool2d(a,2,ceil_mode=True))
+        b = F.interpolate(b,size=a.shape[-2:],mode='bilinear',align_corners=False)
+        return self.head(self.dec(torch.cat((a,b),1)))
+
+class CandidateGate(nn.Module):
+    def __init__(self):
+        super().__init__()
+        self.net=nn.Sequential(nn.Linear(5,16),nn.ReLU(),nn.Linear(16,1))
+    def forward(self,x):
+        return self.net(x).squeeze(-1)
+
+def segmentation_loss(logits, target, name='bce_dice'):
+    p = logits.sigmoid()
+    bce = F.binary_cross_entropy_with_logits(logits,target)
+    dims=tuple(range(1,p.ndim)); eps=1e-6
+    inter=(p*target).sum(dims)
+    if name=='bce': return bce
+    if name=='focal':
+        pt=p*target+(1-p)*(1-target)
+        weight=.25*target+.75*(1-target)
+        return (weight*(1-pt).pow(2)*F.binary_cross_entropy_with_logits(logits,target,reduction='none')).mean()
+    dice=(1-(2*inter+eps)/(p.sum(dims)+target.sum(dims)+eps)).mean()
+    if name=='dice': return dice
+    if name=='bce_dice': return bce+dice
+    if name=='iou': return (1-(inter+eps)/(p.sum(dims)+target.sum(dims)-inter+eps)).mean()
+    raise ValueError(name)
+
+def gate_loss(logits, labels, areas, dual, alpha, size_weight=True, tau=1.):
+    import math
+    pos=labels>0.5; neg=~pos
+    weight=areas.clamp_min(1).rsqrt().clamp(.1,1) if size_weight else torch.ones_like(areas)
+    positive=(F.softplus(-logits[pos])*weight[pos]).sum()/max(1,int(pos.sum()))
+    upper=F.softplus(logits[neg]/tau).sum()/math.log(2)
+    return positive+dual*(upper-alpha), upper