Modern scientific studies routinely collect multiple related responses under common high- dimensional inputs, where each response may depend on predictors through nonlinear effects while sharing latent structure across responses. Existing methods address either nonparamet- ric additive regression for a single response or low-rank regularization in multivariate linear models, but not both simultaneously. We propose a multivariate additive spline framework combining P-spline smoothness regularization with weighted nuclear-norm shrinkage across responses. A central contribution is an efficient two-step fitting procedure: (i) update a data- adaptive reweighting matrix from the current singular structure, and (ii) solve a quadratic surrogate admitting a closed-form ridge-type solution. This formulation enables explicit bias–variance characterization and asymptotic inference, supported by extensive simulations. Two real-data applications, including a CMP removal-rate analysis, further demonstrate the method’s advantages in stability, interpretability, and predictive performance.