{"id":"9287d4fb-bfde-4123-be8d-38c787418d22","arxiv_id":"2507.14631","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":8.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":0,"one_line_summary":"A deterministic polynomial-time algorithm computes a sqrt(d)-approximation to the k-subspace median for every k, avoiding the exponential-in-k runtime of earlier randomized methods.","lead":"A new algorithm finds a robust low-dimensional subspace for point sets, minimizing sum of distances rather than squared distances. It is the first deterministic polynomial-time method whose running time and approximation factor do not grow exponentially with the subspace dimension.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Approximation proof assumes exact SDP minimizer, but the supplied solver is only epsilon-approximate, so the stated sqrt(d) guarantee is not established.","rationale":"I read the paper's central claim as the combination of the sqrt(d) approximation factor and polynomial runtime. The reader's verdict flags Lemma 7.3's false t0 > 1 inequality; that is a real and important bug in the runtime proof. However, I find a more basic gap in the correctness proof itself: Algorithm 1's Line 1 is analyzed with an exact minimizer of the convex relaxation, while Section 7 supplies only an epsilon-approximate solver. Since the additive epsilon error carries through the rounding proof and cannot be eliminated by choosing epsilon small in the zero-optimum case, the theorem as stated does not follow. This concern is independent of the runtime issue and would remain even if Lemma 7.3 were fixed. The algorithm may be salvageable by stating the guarantee as (sqrt(d)+delta)-approximation for any delta > 0, or by proving the rounding is exact on degenerate inputs, but neither is present. In good faith, I do not see a fatal flaw in the approximation factor for nondegenerate inputs; the proof structure is sound if the SDP is solved exactly. Hence the correct verdict remains CONDITIONAL, with the additional condition that the SDP accuracy must be incorporated into the approximation guarantee.","tokens_in":23205,"tokens_out":24090,"duration_ms":284763,"concrete_test":"Analytical check: rerun the proof of Theorem 6.2 with an epsilon-optimal feasible solution X~ in place of the exact minimizer X*; verify that (26) becomes sum_i dist(p_i,S) <= sqrt(d)(relaxksm + epsilon), so for P = {(1,0),(2,0)}, k = 1 (ksm = 0), the certificate would be <= sqrt(2)*epsilon > 0, failing the required <= 0. Computational check: run Algorithm 1 on this input with any SDP solver set to tolerance epsilon = 1e-6; if the returned subspace has positive cost, the stated sqrt(d)-approximation is violated; if it has zero cost, the equality (26) needed by the proof still does not follow from the approximate solver.","verdict_should_be":"CONDITIONAL","load_bearing_attack":"The central guarantee (Corollary 7.4.1) states that Algorithm 1 returns a sqrt(d)-approximation to the k-subspace median in (nd log Delta)^{O(1)} time. The correctness proof of Theorem 6.2 relies on (26), where X* is the exact minimizer of (6) computed in Line 1. However, Section 7 (Proposition 7.2) only establishes that the central-path method computes an additive epsilon-approximation to Problem (6): a feasible (X~, y~) with 1^T y~ <= relaxksm + epsilon. Substituting X~ for X* in the proof yields, instead of (26), sum_i dist(p_i,S) <= sqrt(d) * 1^T y~ <= sqrt(d)(relaxksm + epsilon) <= sqrt(d)*ksm + sqrt(d)*epsilon. For inputs with ksm = 0 (e.g., P = {(1,0),(2,0)}, k = 1), the bound is sqrt(2)*epsilon > 0, so the claimed sqrt(d)-approximation is not certified for any fixed positive epsilon. Corollary 7.4.1 does not specify how epsilon should be chosen or how the additive term is absorbed; in particular, no positive epsilon can satisfy the definition of alpha-approximation when the optimum is zero. Since the central path method converges only asymptotically and no exact SDP solver is provided, the algorithm as stated does not meet its own approximation guarantee. This is a gap in the proof of the main result, and it is more fundamental than the t0 > 1 issue in Lemma 7.3, which only affects the runtime analysis.","agreement_with_reader":"disagree"},"referee_report":{"model":"deepseek-v4-flash","summary":"This paper studies the k-subspace median (kSM) problem, minimizing the sum of Euclidean distances from n points to a k-dimensional linear subspace. The authors relax the non-convex rank-constrained problem (5) to a convex SDP-SOCP problem (6), solve the relaxation, and round the resulting PSD matrix via an eigendecomposition followed by a simple LP-based selection of d-k coordinates (Proposition 6.1). Theorem 6.2 proves, conditionally on computing an exact minimizer of the relaxation, that the rounding yields a sqrt(d)-approximation to kSM. Section 7 proposes a central-path method and claims polynomial running time, leading to Corollary 7.4.1. The paper also includes open-source code and experiments on two UCI datasets.","tokens_in":23501,"tokens_out":11687,"duration_ms":131223,"significance":"If the main result were fully established, it would be a significant advance: a deterministic polynomial-time approximation for kSM with approximation factor sqrt(d) and running time not exponential in k, answering an open problem stated in the paper. The relaxation-rounding idea is elegant, the inequality relaxksm(P,k) <= ksm(P,k) in (20) is correct by feasibility, and Theorem 6.2's proof is structurally sound under the exact-minimizer assumption. The paper also ships open code and reports experiments. However, the gap between the exact minimizer assumed in the correctness proof and the epsilon-approximate solver supplied in Section 7 is load-bearing, and the runtime analysis in Lemma 7.3 contains an invalid inequality; the advertised result is therefore not established in the present form.","major_comments":[{"comment":"The algorithm's own solver does not certify the asserted multiplicative guarantee. Theorem 6.2 uses Eq. (26), which requires that (X*, y*) in Algorithm 1, Line 1 is an exact minimizer of Problem (6). Proposition 7.2 only produces a feasible (X~, y~) with 1^T y~ <= relaxksm(P,k) + epsilon. Repeating the proof with X~ gives sum_i dist(p_i,S) <= sqrt(d) * 1^T y~ <= sqrt(d)(relaxksm(P,k)+epsilon) <= sqrt(d)*ksm(P,k)+sqrt(d)*epsilon. For inputs with ksm(P,k)=0, e.g., P={(1,0),(2,0)} in R^2 with k=1, this bound is sqrt(2)*epsilon > 0, so the output is not certified to be a sqrt(d)-approximation under the paper's own definition in Section 4. Corollary 7.4.1 does not specify how epsilon should be chosen or how the additive term is absorbed, and Theorem 7.4's displayed time bound in (51) omits epsilon altogether. This is a gap in the central claim, not a local typo.","section":"Sec. 7, Algorithm 1, Corollary 7.4.1"},{"comment":"The bound t0 >= 2/sqrt(e) has the wrong direction. Since y0_i = sqrt(||X0 p_i||^2 + e) >= sqrt(e), the denominator of (34) is at least n*sqrt(e), so t0 <= 2/sqrt(e) ~ 1.21, and for large coordinates t0 can be much smaller than 1 (e.g., d=2, k=1, p_1=100 e_1 gives t0 ~ 0.04). The proof then divides by t0-1 in (50), which is invalid when t0 <= 1. Consequently the stated polynomial-time bound on the initial Newton iterations is not established. This does not affect the conditional sqrt(d) approximation proof, but it does affect the running-time half of the main theorem.","section":"Lemma 7.3, Eqs. (34)-(50)"}],"minor_comments":[{"comment":"In (21) the projection onto the output subspace is written as p_i - (I-E)V p_i, but the actual projection matrix defined in Line 10 is V(I-E)V^T. With V orthogonal the intended equality to ||E p_i||_2 becomes correct after replacing (I-E)V by V(I-E)V^T; the surrounding argument should be corrected accordingly.","section":"Theorem 6.2, Eq. (21)"},{"comment":"The displayed chain 'sqrt(d) ||y||_1 <= ||y||_2 <= ||y||_1' is reversed; the correct statement is ||y||_2 <= ||y||_1 <= sqrt(d)||y||_2. The inequalities actually used in (17)-(19) are correct, so this is a presentation error.","section":"Theorem 6.2, after Eq. (19)"},{"comment":"The abstract speaks of approximating the 'affine k-subspace mean', but Problem (5) and the algorithm concern linear subspaces through the origin; the affine case needs the standard reduction mentioned in Section A. Please use the terminology consistently.","section":"Abstract and Section 4"},{"comment":"The comment 'See Proposition 7.2 and Lemma 7.3' next to 'Compute a minimizer' is misleading because Proposition 7.2 only provides an epsilon-approximate solution. The pseudocode and the theorem statements should distinguish clearly between the exact minimizer used in Theorem 6.2 and the approximate solution used in the running-time section.","section":"Algorithm 1, Line 1"}],"recommendation":"major_revision","confidential_remarks":"The paper is a serious attempt and the core rounding scheme is plausible, but the exact-versus-approximate mismatch is fundamental: the current theorem statements describe a guarantee the supplied algorithm does not prove. I would encourage the authors to either provide an exact or properly absorbing treatment of the zero-optimum case, or restate the main result as a bicriteria (sqrt(d), epsilon)-approximation with explicit epsilon-dependence in the running time. The t0 issue in Lemma 7.3 is likely repairable with a different choice of the starting parameter, but it should be fixed before the polynomial-time claim is made."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"The genuinely new thing here is the relaxation: instead of the non-convex rank-(d-k) projection constraint, they minimize the sum of Euclidean distances over the convex hull {0 ⪯ X ⪯ I, trace(X)=d-k}, then round the optimal PSD matrix by keeping the d-k eigenvectors with the largest sum of absolute coordinate projections. That is a clever idea, and the proof that the rounding gives a sqrt(d) approximation — assuming you have an exact minimizer of the relaxation — is short and, as far as I can tell, correct. The feasibility argument for relaxksm ≤ ksm is also clean. The novelty claim relative to the prior exponential-in-k algorithms looks right.\n\nBut the headline theorem is not established. The paper's own Section 7 only shows how to compute an additive epsilon-approximation to the relaxed problem via central path. Corollary 7.4.1 then combines Theorem 6.2 with that epsilon-solver as if Line 1 of Algorithm 1 produced the exact minimizer X*. It does not. Substituting the approximate X~ into the rounding proof gives sqrt(d)*relaxksm + sqrt(d)*epsilon, and when ksm=0 — e.g., points lying on a line, k=1 — the claimed multiplicative sqrt(d) approximation is simply not certified for any fixed positive epsilon. This is not a cosmetic detail; it breaks the main guarantee.\n\nThe second soft spot is in the runtime analysis. Lemma 7.3 assumes t0 > 1 to divide by (t0−1), but t0 = 2n / sum_i sqrt(||X0 p_i||^2 + e) can be as small as O(1/sqrt(e)) for large-coordinate inputs, and the example with p_1 = 100 e_1 in d=2, k=1 gives t0 ≈ 0.04. So the bound on the number of inner Newton iterations is not derived. This only affects the polynomial-time claim, but since that is half the paper's stated contribution, it matters.\n\nThe experiments are a reasonable sanity check but not decisive: they compare against SVD and a single run of a randomized method, and the reported speedups, while real, do not compensate for the proof gaps.\n\nBottom line: this paper deserves a serious referee, but the current version is conditional. The relaxation and rounding technique is worth preserving, and I expect the epsilon issue can be fixed either by adding an exactness assumption, by stating a bicriteria bound, or by a separate zero-optimum case analysis. The t0 issue needs a corrected initial-point argument. Send it to review, but the referees should insist on both fixes before acceptance.","headline":"The relaxation idea is real and the sqrt(d) proof is structurally sound, but the main theorem is not proven as stated because the solver is only epsilon-approximate and the t0 initial-point bound in Lemma 7.3 is false.","tokens_in":24053,"tokens_out":1696,"would_cite":false,"duration_ms":21853,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["68W25","90C22","90C25"],"pacs":[],"model":"deepseek-v4-flash","headline":"The paper proves that the k-subspace median—the subspace minimizing the sum of ordinary, non-squared Euclidean distances—admits a deterministic polynomial-time sqrt(d)-approximation, for every k up to d−1.","keywords":["k-subspace median","mixed norm ℓ_{2,1}","semidefinite programming","second-order cone programming","deterministic approximation algorithm","projection matrix rounding","robust subspace approximation","non-squared Euclidean distances"],"falsifier":"Take d=2, k=1, and a single input point p_1 = 100 e_1. With X_0 = (d−k)/d · I = I/2, compute t_0 = 2n / sum_i $\\sqrt$(||X_0 p_i||^2 + e) = 2 / $\\sqrt$(2500 + e) ≈ 0.04 < 1. Lemma 7.3 divides by t_0−1, so this input directly tests whether the claimed iteration bound covers all allowed inputs; observing the algorithm's termination time on such an input would settle whether Theorem 7.4 holds as stated.","tokens_in":22946,"feed_emoji":"📐","tokens_out":7357,"duration_ms":70744,"temperature":0.7,"pith_summary":"The paper aims to solve the k-subspace median problem: given n points in R^d and an integer k, find the k-dimensional subspace that minimizes the sum of ordinary (non-squared) Euclidean distances to the points. This is a non-convex problem for k < d−1, and previous provable algorithms were randomized with running time or approximation factor exponential in k. The paper claims the first deterministic polynomial-time algorithm for the problem, with approximation factor $\\sqrt$(d) and running time (nd log Δ)^{O(1)} for points with coordinates bounded by Δ. A sympathetic reader should care because the median subspace is more resistant to outliers than PCA, and a polynomial-time guarantee independent of k opens the door to large-k and high-dimensional applications.","feed_headline":"Median subspace: sqrt(d)-approximation in polynomial time","feed_subtitle":"The first deterministic algorithm whose runtime and error do not grow exponentially in k.","key_machinery":"The load-bearing object is the mixed SDP-SOCP relaxation (Problem 6): a convex program whose feasible set is the convex hull of rank-(d−k) projection matrices, combining a positive-semidefinite constraint on X with second-order cone constraints ||X p_i||_2 ≤ y_i. The rounding step uses the vector q ∈ R^d of summed absolute projections onto the eigenaxes of the relaxed solution; choosing the d−k smallest entries of q is exactly a linear-programming min over binary vectors with that many ones (Proposition 6.1). The chain of inequalities ||diag(ζ) V p_i||_2 ≤ ||diag(ζ) V p_i||_1 ≤ ... ≤ $\\sqrt$(d) ||X^* p_i||_2 connects the rounded projection to the relaxation optimum and yields the $\\sqrt$(d) factor.","core_discovery":"The central discovery is that the non-convex k-subspace median can be approximated by a convex relaxation over symmetric matrices: minimize the sum of auxiliary variables y_i subject to 0 ⪯ X ⪯ I, trace(X) = d−k, and ||X p_i||_2 ≤ y_i for every point. The minimizer X* is not itself a projection matrix, but the paper shows how to round it: diagonalize X*, compute q_j = sum_i |(V^T p_i)_j| for each eigenaxis, keep the d−k axes with smallest total absolute projection, and take the resulting projection matrix. This rounding yields a genuine k-subspace whose sum of non-squared distances is at most $\\sqrt$(d) times the optimum, and the relaxation is solvable in polynomial time by a central-path interior-point method. Theorem 6.2 proves the approximation guarantee and Corollary 7.4.1 states the full polynomial running time.","pith_inferences":["Because the approximation factor depends only on d and not on k, the guarantee stays sqrt(d) even when k is Ω(d); for very high-dimensional inputs, composing the method with a dimension-reduction step would trade this factor for one that grows with k instead.","The rounding rule suggests a deterministic selection principle in the eigenbasis of the relaxed solution: pick the axes with smallest total absolute projection. This principle could be tested as a standalone heuristic on noisy datasets without solving the full SDP to high precision.","An immediate empirical check is whether, on datasets with outliers, the returned subspace achieves lower sum-of-distances loss than PCA/SVD, matching the comparison reported in the paper's experiments."],"forward_implications":["For every point set with coordinates in [−Δ,Δ]^d and every k ∈ [1,d−1], the algorithm returns a k-subspace S with sum_i dist(p_i,S) ≤ sqrt(d) · ksm(P,k), in (nd log Δ)^{O(1)} time.","This settles, affirmatively, the open question of whether a deterministic polynomial-time approximation with a factor not exponential in k exists for the k-subspace median.","The same algorithm can be run on coresets or sketches of the input to obtain near-linear time in n, at the price of randomization, while keeping the approximation factor.","The techniques are expected to extend to other mixed ℓ_{2,z} objectives, such as the k-subspace center (z=∞), and to settings with outliers and sparsity demands."],"supporting_citations":[{"why":"Supplies the main baseline: the prior randomized (1+ε)-approximation whose running time and success probability are exponential in k.","marker":"[17]"},{"why":"Provides the sampling-based dimension reduction that produces a small point subset whose span contains a (1+ε)-approximation, a key comparison point for the new approach.","marker":"[19]"},{"why":"Gives input-sparsity algorithms and NP-hardness for robust subspace approximation, framing the polynomial-time question and the benchmark.","marker":"[21]"},{"why":"The low-rank approximation theorem that the rounding step adapts from Frobenius error to the mixed ℓ_{2,1} norm.","marker":"[29]"},{"why":"The Courant–Fischer selection principle invoked alongside the low-rank approximation as the classical counterpart of the new rounding rule.","marker":"[30]"},{"why":"A convex-optimization reference that supplies the PSD and second-order cones, self-concordant barriers, and interior-point background used in the central path analysis.","marker":"[31]"},{"why":"Supplies the central path method and its iteration bounds used to solve the convex relaxation in polynomial time.","marker":"[34]"}],"fun_headline_variants":["Non-squared Euclidean k-subspace median: sqrt(d)-approx in poly time","First poly-time sqrt(d)-approx for k-subspace median","Rounding a convex relaxation achieves sqrt(d)-approx median subspace","Sparse robust k-subspace median: polynomial algorithm, sqrt(d) error","Deterministic poly-time algorithm for non-squared distances k-subspace"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The polynomial-time part of the proof assumes that the initial-point parameter t0 defined in Lemma 7.3 is always larger than 1, because the argument divides by t0−1; for inputs with large coordinates t0 can be smaller than 1, so the stated bound on the number of interior-point iterations is not derived.","fun_headline_variants_meta":{"raw":{"variants":["Non-squared Euclidean k-subspace median: sqrt(d)-approx in poly time","First poly-time sqrt(d)-approx for k-subspace median","Rounding a convex relaxation achieves sqrt(d)-approx median subspace","Sparse robust k-subspace median: polynomial algorithm, sqrt(d) error","Deterministic poly-time algorithm for non-squared distances k-subspace"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000205,"raw_usage":{"total_tokens":1431,"prompt_tokens":1018,"completion_tokens":413,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":634,"completion_tokens_details":{"reasoning_tokens":321}},"tokens_in":634,"tokens_out":413,"duration_ms":4564,"temperature":1.0,"reasoning_tokens":321,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-06T15:53:07.944207+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Take d=2, k=1, and a single input point p_1 = 100 e_1. With X_0 = (d−k)/d · I = I/2, compute t_0 = 2n / sum_i $\\sqrt$(||X_0 p_i||^2 + e) = 2 / $\\sqrt$(2500 + e) ≈ 0.04 < 1. Lemma 7.3 divides by t_0−1, so this input directly tests whether the claimed iteration bound covers all allowed inputs; observing the algorithm's termination time on such an input would settle whether Theorem 7.4 holds as stated.","supporting_citations":[{"cited_title":"Efficient subspace approximation algorithms,","cited_arxiv_id":null,"evidence_quote":"Supplies the main baseline: the prior randomized (1+ε)-approximation whose running time and success probability are exponential in k."},{"cited_title":"Sampling-based dimension reduction for subspace approximation,","cited_arxiv_id":null,"evidence_quote":"Provides the sampling-based dimension reduction that produces a small point subset whose span contains a (1+ε)-approximation, a key comparison point for the new approach."},{"cited_title":"Input sparsity and hardness for robust subspace approximation,","cited_arxiv_id":null,"evidence_quote":"Gives input-sparsity algorithms and NP-hardness for robust subspace approximation, framing the polynomial-time question and the benchmark."},{"cited_title":"On least squares solutions subject to a rank restriction,","cited_arxiv_id":null,"evidence_quote":"The low-rank approximation theorem that the rounding step adapts from Frobenius error to the mixed ℓ_{2,1} norm."},{"cited_title":"The monotonicity theorem, cauchy’s interlace theorem, and the courant-fischer theorem,","cited_arxiv_id":null,"evidence_quote":"The Courant–Fischer selection principle invoked alongside the low-rank approximation as the classical counterpart of the new rounding rule."},{"cited_title":"Interior point polynomial methods in convex programming,","cited_arxiv_id":null,"evidence_quote":"Supplies the central path method and its iteration bounds used to solve the convex relaxation in polynomial time."}],"review_version":1}