Forum for Interdisciplinary Mathematics

Xiang-Sheng Wang
Chisheng Wang

Machine
Learning
in Data
Processing
Forum for Interdisciplinary Mathematics

Editor-in-Chief
Viswanath Ramakrishna, University of Texas, Richardson, USA

Editorial Board Members

Ashis SenGupta, Indian Statistical Institute, Kolkata, India
Balasubramaniam Jayaram, Indian Institute of Technology, Hyderabad, India
P. V. Subrahmanyam, Indian Institute of Technology Madras, Chennai, India
Ravindra B. Bapat, Indian Statistical Institute, New Delhi, India
The Forum for Interdisciplinary Mathematics is a Scopus-indexed book series.
It publishes high-quality textbooks, monographs, contributed volumes and lecture
notes in mathematics and interdisciplinary areas where mathematics plays a fun-
damental role, such as statistics, operations research, computer science, financial
mathematics, industrial mathematics, and bio-mathematics. It reflects the increasing
demand of researchers working at the interface between mathematics and other
scientific disciplines.
Xiang-Sheng Wang • Chisheng Wang

Machine Learning in Data

Processing
Xiang-Sheng Wang Chisheng Wang
University of Louisiana at Lafayette Shenzhen University
Lafayette, LA, USA Shenzhen, Guangdong, China

ISSN 2364-6748 ISSN 2364-6756 (electronic)

Forum for Interdisciplinary Mathematics
ISBN 978-3-032-20854-5 ISBN 978-3-032-20855-2 (eBook)
[Link]

Mathematics Subject Classification: 62J07, 15A18, 65K05, 68T07, 90C25

© The Editor(s) (if applicable) and The Author(s), under exclusive license to Springer Nature Switzerland
AG 2026

This work is subject to copyright. All rights are solely and exclusively licensed by the Publisher, whether
the whole or part of the material is concerned, specifically the rights of translation, reprinting, reuse
of illustrations, recitation, broadcasting, reproduction on microfilms or in any other physical way, and
transmission or information storage and retrieval, electronic adaptation, computer software, or by similar
or dissimilar methodology now known or hereafter developed.
The use of general descriptive names, registered names, trademarks, service marks, etc. in this publication
does not imply, even in the absence of a specific statement, that such names are exempt from the relevant
protective laws and regulations and therefore free for general use.
The publisher, the authors and the editors are safe to assume that the advice and information in this book
are believed to be true and accurate at the date of publication. Neither the publisher nor the authors or
the editors give a warranty, expressed or implied, with respect to the material contained herein or for any
errors or omissions that may have been made. The publisher remains neutral with regard to jurisdictional
claims in published maps and institutional affiliations.

This Springer imprint is published by the registered company Springer Nature Switzerland AG
The registered company address is: Gewerbestrasse 11, 6330 Cham, Switzerland

If disposing of this product, please recycle the paper.

To Dongying Xu
Preface

Machine learning is perceived in diverse ways across academic and professional

communities, with each perspective influenced by disciplinary background and
practical experience. Rooted in mathematics, statistics, computer science, data
science, and numerical computation, the field has advanced rapidly in the era of
big data. Its broad impact on science, engineering, and everyday life has led many
universities to actively integrate machine learning into both undergraduate and
graduate curricula.
At the request of several graduate students, the first author began teaching a
course on machine learning in Fall 2022 at the University of Louisiana at Lafayette.
Offered through the Department of Mathematics, the course attracted students from
a broad range of disciplines, including mechanical engineering, systems engineer-
ing, electrical engineering, earth and energy sciences, and industrial chemistry.
Many of these students had substantial hands-on experience applying machine
learning techniques in data processing but wished to gain a deeper understanding
of the underlying mathematical principles that explain why these methods work so
effectively.
In contrast, most mathematics majors enrolled in the course had little or no
prior exposure to machine learning or computer programming and often lacked
experience writing even simple code. Their interest was motivated by the growing
prominence of machine learning and a desire to understand foundational concepts
such as neural networks and support vector machines from a mathematical stand-
point.
To accommodate this diverse student population, the first author consulted with
the second author, a geo-informatics researcher specializing in the application of
mathematical and data analysis techniques to geoscience problems. Together, the
authors reviewed numerous textbooks and online resources. While many of these
materials are well suited for engineering students focused on applications, they
are often too advanced for mathematics students with limited programming back-
ground. Moreover, existing resources tend to emphasize large real-world datasets
and extensive coding, often at the expense of clear and systematic mathematical
exposition.

vii
viii Preface

This lack of suitable instructional material motivated the authors to write a

textbook specifically designed for students with minimal programming experience.
The goal of the book is to provide a mathematically grounded introduction to
machine learning without requiring extensive coding skills.
During the 2023–2024 academic year, the first author took a sabbatical and
collaborated closely with the second author to substantially revise and expand
the lecture notes from the Fall 2022 course. This collaboration resulted in the
completion of the initial draft of the book. Upon returning from sabbatical, the first
author used the manuscript as the primary text for a redesigned machine learning
course in Fall 2024. The chapters were structured to fit a standard one-semester
course of 30–40 contact hours, and an appendix was added to review essential topics
in calculus and linear algebra.
This book consolidates fundamental mathematical derivations and theoretical
results that are otherwise scattered throughout the machine learning literature.
Unlike many existing texts that focus on large-scale real-world datasets, this book
emphasizes simplified synthetic data to clearly illustrate core concepts and methods.
While coding is optional, concise Python scripts are provided to allow students to
experiment with examples and modify them for exercises.
The only prerequisites for this book are a solid understanding of calculus and
linear algebra. Key topics such as Taylor expansions, derivatives, the chain rule, and
matrix operations are used extensively. For example, the derivation of backward
propagation in neural networks relies heavily on the chain rule and matrix calculus.
Basic statistical concepts—such as mean, variance, probability distributions, and
probability density functions—are occasionally referenced, and familiarity with
elementary statistics will aid comprehension.
This book focuses primarily on supervised learning, where models are trained
using labeled data, consisting of input–output pairs, and their performance is
evaluated by their ability to generalize from training data to previously unseen
test data. While this book does not treat them in detail, we emphasize that other
major paradigms of machine learning also exist, including unsupervised learning,
which seeks structure in unlabeled data, and reinforcement learning, which studies
decision-making through interaction with an environment.
Throughout the book, learning is formulated as an optimization problem, and
model performance must be assessed in terms of prediction accuracy and generaliza-
tion error. Issues such as overfitting are discussed from a mathematical perspective,
motivating regularization techniques and model selection. The text is intentionally
selective: it concentrates on regression models, classification methods, support
vector machines, and neural networks, while omitting topics such as probabilistic
graphical models and reinforcement learning in order to maintain a coherent and
rigorous mathematical narrative. The goal is not to provide a comprehensive survey
of all machine learning methods, but to equip readers—particularly those with
strong mathematical or computational background—with a clear understanding of
the numerical methods and theoretical principles that underlie modern supervised
learning.
Preface ix

This book is intended for two primary audiences: students majoring in math-
ematics or statistics who are not actively working in machine learning but wish
to understand its fundamental concepts and terminology, and students from other
disciplines who have practical experience using machine learning tools and seek
a rigorous understanding of the mathematical foundations underlying these tech-
niques.
A central objective of this book is to promote machine learning as a standard
undergraduate or graduate course—comparable to differential equations—requiring
minimal prerequisites in programming and statistics. By relying primarily on cal-
culus and linear algebra, the material is accessible to students who have completed
basic coursework in these subjects.
To achieve this accessibility, the presentation is deliberately elementary and self-
contained. Upon completing this book or a course based on it, students should
acquire foundational knowledge of linear and nonlinear regression, regularization,
neural networks, batch normalization, support vector machines, gradient-based
optimization methods, and principal component analysis. The exposition follows
a theorem-and-proof style to ensure mathematical rigor and clarity.
Python is chosen as the programming language because it is free, widely
available, and easy to read. Although programming languages evolve over time,
Python’s simplicity makes it particularly suitable for beginners. Code examples are
intentionally kept short—typically no more than 30–40 lines—and avoid reliance
on specialized libraries.
The first author gratefully acknowledges the unwavering support of his wife,
Hongying, and his children, Mingqian and Tongwei. The second author gratefully
acknowledges the constant encouragement of his wife, Yanhong, and his children,
Yuchen and Xihe. The authors are deeply thankful for their families’ support, which
enabled them to persevere through the challenges of completing this book.

Lafayette, LA, USA Xiang-Sheng Wang

Shenzhen, China Chisheng Wang
Contents

1 Matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.1 Notations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2 Determinant. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3 Block Matrix. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.4 Rank. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.5 Fréchet Derivative . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2 Linear Regression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.1 Notations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.2 Linear Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.3 Maximum Likelihood Estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.4 Least Squares Approximation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.5 Sum of Squared Errors . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.6 Variance Inflation Factor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.7 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3 Regularization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.1 Overfitting Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.2 Ridge Regression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.3 Convex Optimization. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.4 LASSO. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.5 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
3.6 Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
4 Nonlinear Regression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4.1 Nonlinear Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4.2 Sigmoid Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
4.3 Optimization Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.4 Gradient Descent Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

xi
xii Contents

4.5 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50

4.6 Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
5 Shallow Neural Network . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
5.1 Notations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
5.2 The XOR Gate. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
5.3 Composite Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
5.4 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
5.5 Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
6 Deep Neural Network . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
6.1 The Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
6.2 Backward Propagation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
6.3 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
6.4 Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
7 Batch Normalization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
7.1 The Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
7.2 Backward Propagation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
7.3 Python Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
7.4 Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
8 Support Vector Machine . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
8.1 Notations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
8.2 Optimization Problems. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
8.3 Equivalence Theorem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
9 Gradient Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
9.1 Conjugate Gradient Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
9.2 One-Step Method. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86
9.3 Two-Step Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
10 Dimensionality Reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93
10.1 Schur Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93
10.2 Singular Value Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
10.3 Principal Component Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100

A Related Topics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103

A.1 Calculus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
A.2 Linear Algebra. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106
Contents xiii

B Hints to Selected Exercise Problems. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111

C Further Readings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113
References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115
Index . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
Chapter 1
Matrix

Abstract In this chapter, we review some terminologies and properties of matrices

which will be used throughout this book. For simplicity, we mainly focus on
real matrices, even though many statements and proofs can be easily extended to
matrices in complex number system.

1.1 Notations

Throughout this chapter, we denote by Rn . the space of n-dimensional real vectors.

By default, a vector v ∈ Rn . is a column vector with components v1 , · · · , vn .. A
collection of m vectors in Rn . forms a matrix A = (A1 , · · · , Am ). of dimension
n × m., where
⎛⎞
A1k
⎜ . ⎟
.Ak = ⎝ . ⎠ ∈ R , k = 1, · · · , m.
n
.
Ank

The space of all n × m. matrices is denoted by Rn × · · · × Rn = Rn×m . (m copies).

Given A ∈ Rn×m . with n > 1. and m > 1., the (i,j)-minor of A is the submatrix
A ) ∈ R(n−1)×(m−1) . obtained by deleting the i-th row and j -th column of A.
(i,j

Explicitly,

(i,j )
Akl
. = Ak l , 1 ≤ k ≤ n − 1, 1 ≤ l ≤ m − 1,

where

k, k < i, l, l < j,
k =
. l =
k + 1, k ≥ i, l + 1, l ≥ j.

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 1

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
2 1 Matrix

The transpose of a matrix A ∈ Rn×m . is the matrix AT ∈ Rm×n . defined by

ATjk = Akj , 1 ≤ j ≤ m, 1 ≤ k ≤ n.
.

In particular, a column vector x ∈ Rn = Rn×1 . can be viewed as an n × 1.

matrix and its transpose x T ∈ R1×n . is a row vector. A square matrix A ∈ Rn×n . is
symmetric if AT = A.; namely, Aj k = Akj . for any 1 ≤ j, k ≤ n.. The trace of a
square matrix A ∈ Rn×n . is defined as
n
. tr(A) = Akk .
k=1

The product of two matrices A ∈ Rn×m . and B ∈ Rm×p . is the matrix AB ∈ Rn×p .
given by
m
(AB)ij =
. Aik Bkj , 1 ≤ i ≤ n, 1 ≤ j ≤ p.
k=1

The standard basis of Rn . consists of vectors e1 , · · · , en ., where the j -th

component of ek . is

1, j = k,
δj k =
.
0, j k.

The identity matrix is I = (e1 , · · · , en ) ∈ Rn×n .. It is easily seen that AI = I A =

A. for any A ∈ Rn×n .. A square matrix A ∈ Rn×m . is invertible (or nonsingular) if
there exists B = A−1 ∈ Rn×m . such that

AB = BA = I.
.

Let (i1 , · · · , in ). be any permutation of the indices (1, · · · , n).. Denote by

σ (i1 , · · · , in ). the minimum number of adjacent transpositions required to convert
(i1 , · · · , in ). into (1, · · · , n).. The sign of the permutation is defined as

s(i1 , · · · , in ) = (−1)σ (i1 ,··· ,in ) .

.

Clearly, the sign function is antisymmetric, meaning that

s(i1 , · · · , ij , · · · , ik , · · · , in ) = −s(i1 , · · · , ik , · · · , ij , · · · , in ).
.
1.2 Determinant 3

1.2 Determinant

Given a function f on Rn×n ., we say that f is multilinear if, for any k = 1, · · · , n.,
any scalars α, β ∈ R., and any vectors A1 , · · · , An , B1 , · · · , Bn ∈ Rn .,

f (A1 , · · · , Ak−1 , αAk + βBk , Ak+1 , · · · , An )

.

=αf (A1 , · · · , Ak−1 , Ak , Ak+1 , · · · , An )

+ βf (A1 , · · · , Ak−1 , Ak , Ak+1 , · · · , An ). (1.1)

We say that f is antisymmetric if

f (A1 , · · · , Aj , · · · , Ak , · · · , An ) = −f (A1 , · · · , Ak , · · · , Aj , · · · , An ),
.

(1.2)
that is, f changes sign whenever two of its vector arguments are interchanged. We
say that f is normalized if

f (I ) = f (e1 , · · · , en ) = 1,
. (1.3)

where e1 , · · · , en . are the standard basis vectors of Rn . and I = (e1 , · · · , en ). is the

identity matrix in Rn×n .. The determinant of a square matrix A ∈ Rn×n . is the unique
function on Rn×n . that is normalized, antisymmetric, and multilinear. We will show
that such a function exists and is unique.
To establish existence, we define f recursively as follows: f (a) = a . for a ∈ R.
and

A, A ∈ R,
f (A) =
.
n
(1.4)
i=1 (−1)
i−1 A f (A(i,1) ),
i1 A ∈ Rn×n , n > 1.

When n = 2. and A ∈ R2×2 ., the formula (1.4) yields

. f (A) = A11 A22 − A21 A12 .

It is straightforward to verify that f defined in this way is normalized, antisymmet-

ric, and multilinear on R2×2 ..
Lemma 1.1 The function f defined recursively by (1.4) is normalized, antisymmet-
ric, and multilinear on Rn×n ..
Proof We proceed by induction on n. For n = 1. and n = 2., it is evident that f is
normalized, antisymmetric, and multilinear.
Assume that f is normalized, antisymmetric, and multilinear on R(n−1)×(n−1) ..
We will show that these properties hold on Rn×n .. Let k = 2, · · · , n., α, β ∈ R., and
A1 , · · · , An , B1 , · · · , Bn ∈ Rn .. It follows from the recursive definition (1.4) that
4 1 Matrix

f (A1 , · · · , Ak−1 , αAk + βBk , Ak+1 , · · · , An )

.

n
= (−1)i−1 Ai1 f (A(i) (i) (i) (i) (i)
2 , · · · , Ak−1 , αAk + βBk , Ak+1 , · · · , An ),
(i)

i=1

where A(i)j ∈ R
n−1 . (respectively, B (i) ∈ Rn−1 .) denotes the vector obtained from
j
Aj . (respectively, Bj .) by deleting its i-th component. Since f is multilinear on
R(n−1)×(n−1) ., we have

f (A1 , · · · , Ak−1 , αAk + βBk , Ak+1 , · · · , An )

.

n
=α (−1)i−1 Ai1 f (A(i) (i) (i) (i)
2 , · · · , Ak−1 , Ak , Ak+1 , · · · , An )
(i)

i=1
n
(i) (i) (i) (i)
+β (−1)i−1 Ai1 f (A2 , · · · , Ak−1 , Bk , Ak+1 , · · · , A(i)
n )
i=1

=αf (A1 , · · · , Ak−1 , Ak , Ak+1 , · · · , An )

+ βf (A1 , · · · , Ak−1 , Bk , Ak+1 , · · · , An ).

Next, we obtain

f (αA1 + βB1 , A2 , · · · , An )
.

n
(i)
= (−1)i−1 (αAi1 + βBi1 )f (A2 , · · · , A(i)
n )
i=1

=αf (A1 , A2 , · · · , An ) + βf (B1 , A2 , · · · , An ).

Coupling the above two results shows that f is multilinear on Rn×n ..

Recall that e1 , · · · , en . are the standard basis vectors of Rn .. It is easily seen that
e2 , · · · , en(1) . form the standard basis vectors of Rn−1 ., where ej(1) . is the vector
(1)

obtained by deleting the first component of ej . for j = 2, · · · , n.. Consequently,

we have
(1)
f (e1 , · · · , en ) = f (e2 , · · · , en(1) ) = 1,
.

which proves that f is normalized on Rn×n ..

Finally, we assume Aj = Ak . for some 1 ≤ j < k ≤ n.. If j > 1., then the
antisymmetry of f on R(n−1)×(n−1) . together with the recursive definition (1.4) gives

f (A1 , · · · , Aj , · · · , Ak , · · · , An ) = 0.
.
1.2 Determinant 5

If j = 1., then by antisymmetry of f with respect to A2 , · · · , An ., we may assume

without loss of generality that k = 2. and A1 = A2 .. For each l = 3, · · · , n.,
(i)(j )
we denote by Al ∈ Rn−2 . the vector obtained by deleting the i-th and j -th
components of Al .. Applying (1.4) twice yields
n n
(i)(j ) (i)(j )
f (A1 , A2 , · · · , An ) =
. (Ai1 Aj 2 − Aj 1 Ai2 )f (A3 , · · · , An ) = 0,
i=1 j =1

since A1 = A2 .. Hence, f is antisymmetric in all its arguments A1 , A2 , · · · , An ..

This completes the proof.
The following lemma establishes the uniqueness of the determinant function.
Lemma 1.2 If f is a normalized, antisymmetric, and multilinear function on Rn×n .,
then we have
n
f (A1 , · · · , An ) =
. s(i1 , · · · , in ) Aij j , (1.5)
(i1 ,··· ,in )∈Pn j =1

where Pn . denotes the set of all permutations of the indices (1, 2, · · · , n)..
Proof By the multilinearity of f , we obtain
n
f (A1 , · · · , An ) =
. f (ei1 , · · · , ein ) Aij j .
(i1 ,··· ,in )∈Pn j =1

Since f is normalized and antisymmetric, it follows that

f (ei1 , · · · , ein ) = s(i1 , · · · , in ),

.

where s(i1 , · · · , in ). denotes the sign of the permutation (i1 , · · · , in ).. This completes
the proof.
According to Lemma (1.1) and Lemma (1.2), the determinant function is
uniquely determined by the conditions of normalization, antisymmetry, and mul-
tilinearity. From now on, we shall use det. to denote the determinant function. The
following properties are easy to verify.
Proposition 1.1 Given A = (A1 , · · · , An ) ∈ Rn×n ., we have

. det(A) = det(AT ). (1.6)

For any i, j = 1, · · · , n., we have

6 1 Matrix

n
det(A), i = j,
. (−1)i+k Aki det(A(k,j ) ) = det(A)δij = (1.7)
k=1 0, i j.

For any α ∈ R., we have

. det(A1 , · · · , Aj , · · · , Ak + αAj , · · · , An ) = det(A1 , · · · , An ), . (1.8)

det(A1 , · · · , Ak−1 , αAk , Ak+1 , · · · , An ) =α det(A1 , · · · , An ). (1.9)

For any A, B ∈ Rn×n ., we have

. det(AB) = det(A) det(B). (1.10)

In particular, if A is invertible, then det(A) 0.. On the other hand, if det(A) 0.,
then A is invertible with

(−1)i+j det(A(j,i) )
(A−1 )ij =
. . (1.11)
det(A)

1.3 Block Matrix

If A ∈ Rm×n ., B ∈ Rm×q ., C ∈ Rp×n ., and D ∈ Rp×q ., we define the block matrix

A B
M=
. ∈ R(m+p)×(n+q)
C D

such that
⎧
⎪
⎪Aij , i = 1, · · · , m, j = 1, · · · , n,
⎪
⎪
⎨B i = 1, · · · , m, j = n + 1, · · · , n + q,
i,j −n ,
.Mij =
⎪
⎪ Ci−m,j , i = m + 1, · · · , m + p, j = 1, · · · , n,
⎪
⎪
⎩
Di−m,j −n , i = m + 1, · · · , m + p, j = n + 1, · · · , n + q.

The operations of addition, subtraction, and scalar multiplication on block

matrices can be carried out blockwise. The following lemma presents the formula
for block matrix multiplication.
Lemma 1.3 If A ∈ Rm×n ., B ∈ Rm×q ., C ∈ Rp×n ., D ∈ Rp×q ., Ã ∈ Rn×r ., B̃ ∈
Rn×s ., C̃ ∈ Rq×r ., and D̃ ∈ Rq×s ., then

A B Ã B̃ AÃ + B C̃ AB̃ + B D̃
. = . (1.12)
C D C̃ D̃ C Ã + D C̃ C B̃ + C D̃
1.3 Block Matrix 7

Proof Let

A B Ã B̃
M=
. ∈ R(m+p)×(n+q) , M̃ = ∈ R(n+q)×(r+s) .
C D C̃ D̃

Then
n+q
(M M̃)ij =
. Mik M̃kj , i = 1, · · · , m + p, j = 1, · · · , r + s.
k=1

For i = 1, · · · , m. and j = 1, · · · , r ., we have

n n+q
(M M̃)ij =
. Aik Ãkj + Bi,k−n C̃k−n,j = (AÃ + B C̃)ij .
k=1 k=n+1

We check each block explicitly. For i = 1, · · · , m. and j = r + 1, · · · , r + s ., we

have
n n+q
(M M̃)ij =
. Aik B̃k,j −r + Bi,k−n D̃k−n,j −r = (AB̃ + B D̃)i,j −r .
k=1 k=n+1

For i = m + 1, · · · , m + p . and j = 1, · · · , r ., we have

n n+q
(M M̃)ij =
. Ci−m,k Ãkj + Di−m,k−n C̃k−n,j = (C Ã + D C̃)i−m,j .
k=1 k=n+1

For i = m + 1, · · · , m + p . and j = r + 1, · · · , r + s ., we have

n n+q
(M M̃)ij =
. Ci−m,k B̃k,j −r + Di−m,k−n D̃k−n,j −r = (C B̃ + D D̃)i−m,j −r .
k=1 k=n+1

This completes the proof.

The following result is useful for computing the determinant of block matrices.
Lemma 1.4 Let

A B
M=
. ∈ R(m+n)×(m+n) ,
C D

where A ∈ Rm×m ., B = 0 ∈ Rm×n ., C ∈ Rn×m ., and D ∈ Rn×n .. Then

8 1 Matrix

A 0
. det(M) = det = det(A) det(D). (1.13)
C D

Proof We prove by induction on m. For m = 1., the result is immediate. Assume

the formula holds for m − 1.. For each j = 1, · · · , m., we denote by M (1,j ) ∈
R(m+n−1)×(m+n−1) . the submatrix obtained by deleting the first row and the j -th
column of M:

A(1,j ) 0
M (1,j ) =
. ,
C (0,j ) D

where A(1,j ) ∈ R(m−1)×(m−1) . is obtained from A by removing its first row and j -th
column, and C (0,j ) ∈ Rn×(m−1) . is obtained from C by removing its j -th column.
By the induction hypothesis,

. det(M (1,j ) ) = det(A(1,j ) ) det(D).

Expanding det(M). along its first row gives

m
. det(M) = (−1)j +1 A1j det(M (1,j ) )
j =1
m
= (−1)j +1 A1j det(A(1,j ) ) det(D)
j =1

= det(A) det(D),

which completes the proof.

Proposition 1.2 Let

A B
M=
. ∈ R2n×2n ,
C D

where A, B, C, D ∈ Rn×n . and AC = CA.. Then

A B
. det(M) = det = det(AD − CB). (1.14)
C D

Proof We first consider the case where A is invertible; that is, det(A) 0.. Observe
that

I 0 A B A B
. = .
−CA−1 I C D 0 D − CA−1 B
1.4 Rank 9

Taking determinants on both sides and applying Proposition 1.1 and Lemma 1.4
give

A B
. det = det(A) det(D − CA−1 B) = det(AD − ACA−1 B).
C D

Since AC = CA., it follows that det(M) = det(AD − CB)..

Next, assume that det(A) = 0.. Consider the polynomial

.p(λ) = det(A + λI ) = λn + · · · , λ ∈ R.

By the fundamental theorem of algebra, p(λ). has at most n real roots. Hence,
there exists ε0 > 0. such that p(ε) 0. for all ε ∈ (0, ε0 ).. For each such ε., set
Aε = A + εI ., which is invertible. By the previous argument, we have

Aε B
. det = det(Aε D − CB).
C D

Letting ε → 0+ . on both sides yields (1.14). This completes the proof.

1.4 Rank

We say that the rank of a matrix A ∈ Rm×n . is r if there exist invertible matrices
P ∈ Rm×m . and Q ∈ Rn×n . such that

P AQ = (e1 , · · · , er , 0, · · · , 0) ∈ Rm×n ,
.

where e1 , · · · , em . are the standard basis vectors in Rm .. The following lemma

ensures that the rank is well defined.
Lemma 1.5 If

. A = P (e1 , · · · , er , 0, · · · , 0) = (e1 , · · · , es , 0, · · · , 0)Q ∈ Rm×n ,

where P ∈ Rm×m . and Q ∈ Rn×n . are invertible, then r = s ..

Proof For any i = 1, · · · , m., we have

Pij , j = 1, · · · , r,
Aij =
.
0, j = r + 1, · · · , n.

Similarly, for any j = 1, · · · , n., we have

10 1 Matrix

Qij , i = 1, · · · , s,
. Aij =
0, i = s + 1, · · · , m.

Combining these two formulas gives

Pij =0, i = s + 1, · · · , m, j = 1, · · · , r,
.

Qij =0, i = 1, · · · , s, j = r + 1, · · · , n.

Suppose r < s .. We can write Q as a block matrix:

B 0
Q=
. ,
C D

where B = (B1 , · · · , Br , 0, · · · , 0) ∈ Rs×s ., C ∈ R(n−s)×s ., and D ∈ R(n−s)×(n−s) ..

By Lemma 1.4, we have det(Q) = det(B) det(D) = 0., contradicting the
invertibility of Q.
Suppose r > s .. We can write P T . as a block matrix:

E 0
.PT = ,
F G

where E = (E1 , · · · , Es , 0, · · · , 0) ∈ Rr×r ., F ∈ R(m−r)×r ., and G ∈

R(m−r)×(m−r) .. By Proposition 1.1 and Lemma 1.4, we have det(P ) = det(P T ) =
det(E) det(G) = 0., which contradicts the invertibility of P .
Hence r = s .. This completes the proof.
The following results are immediate consequences of the definition.
Proposition 1.3 Let A ∈ Rm×n .. Then

. rank(A) = rank(AT ) ≤ min{m, n}.

Moreover, for any A, B ∈ Rm×n ., we have rank(A) = rank(B). if and only if

there exist invertible matrices P ∈ Rm×m . and Q ∈ Rn×n . such that A = P BQ..
In particular, a square matrix A ∈ Rn×n . is invertible if it has full rank, that is,
rank(A) = n..
We are now ready to show that matrix multiplication cannot increase the rank.
Proposition 1.4 For any A ∈ Rn×m . and B ∈ Rm×p ., we have

. rank(AB) ≤ min{rank(A), rank(B)}. (1.15)

Proof For any matrix M = P (e1 , · · · , er , 0, · · · , 0)Q., where P and Q are

invertible, we have
1.4 Rank 11

Q 0
(M, 0) = P (e1 , · · · , er , 0, · · · , 0)
. ,
0 I

and hence

. rank(M, 0) = rank(M).

For any R = (R1 , · · · , Rm ) ∈ Rn×m . and any s ≤ m., we have

. rank(R(e1 , · · · , es , 0, · · · , 0)) = rank(R1 , · · · , Rs , 0, · · · , 0)

= rank(R1 , · · · , Rs ) ≤ s.

Now assume that B = P (e1 , · · · , es , 0, · · · , 0)Q., where P and Q are invertible.

Then

. rank(AB) = rank(AP (e1 , · · · , es , 0, · · · , 0)) ≤ s = rank(B).

Furthermore

. rank(AB) = rank((AB)T ) = rank(B T AT ) ≤ rank(AT ) = rank(A).

This completes the proof.

The vectors A1 , · · · , Ak ∈ Rn . are said to be linearly dependent if there exists a
nonzero vector v ∈ Rk . such that

(A1 , · · · , Ak )v = A1 v1 + · · · + Ak vk = 0.
.

They are said to be linearly independent if the vector equation

(A1 , · · · , Ak )v = 0
.

has only the trivial solution v = 0 ∈ Rk ..

Lemma 1.6 The vectors A1 , · · · , Ak ∈ Rn . are linearly independent if and only if
rank(A1 , · · · , Ak ) = k..
Proof Suppose (A1 , · · · , Ak )v = 0. for some nonzero vector v ∈ Rk .. Then there
exists an index i ∈ {1, · · · , k}. such that vi 0.. Define

Q = (e1 , · · · , ei−1 , ei+1 , · · · , ek , v) ∈ Rk×k ,

.

where e1 , · · · , ek . denote the standard basis vectors of Rk .. Since det(Q) = vi 0.,

the matrix Q is invertible. Note that
12 1 Matrix

(A1 , · · · , Ak )Q = (A1 , · · · , Ai−1 , Ai+1 , · · · , Ak , 0),

.

which implies

. rank(A1 , · · · , Ak ) = rank(A1 , · · · , Ai−1 , Ai+1 , · · · , Ak ) ≤ k − 1 < k.

Conversely, if rank(A1 , · · · , Ak ) = r < k ., then there exist invertible matrices

P = (P1 , · · · , Pn ) ∈ Rn×n . and Q = (Q1 , · · · , Qk ) ∈ Rk×k . such that

(A1 , · · · , Ak )Q = P (e1 , · · · , er , 0, · · · , 0) = (P1 , · · · , Pr , 0, · · · , 0).

.

Let v = Qr+1 ∈ Rk .. Then (A1 , · · · , Ak )v = 0., where v 0. since det(Q) 0..

Thus, the vectors A1 , · · · , Ak ∈ Rn . are linearly dependent. This completes the
proof.
Let A = (A1 , · · · , Am ) ∈ Rn×m .. Then rank(A) = r . if and only if there exist
linearly independent vectors P1 , · · · , Pr ∈ Rn . such that each column A1 , · · · , Am .
of A can be expressed as a linear combination of P1 , · · · , Pr ..

1.5 Fréchet Derivative

A function f ∈ C(Rn , R). is said to be Fréchet differentiable at x ∈ Rn . if there

exists a bounded linear operator L ∈ C(Rn , R). such that

|f (x + h) − f (x) − Lh|
. lim = 0.
h→0 h

Any bounded linear operator L ∈ C(Rn , R). can be represented by a column

vector v ∈ Rn . satisfying Lh = v T h. for all h ∈ Rn .. In this case, we denote f (x) =
v ∈ Rn . and call f (x). the gradient of f at x ∈ Rn .. It follows that
⎛ ∂f ⎞
∂x1 (x)
⎜ ⎟
f (x) = ⎝ · · · ⎠ ∈ Rn .
. (1.16)
∂f
∂xn (x)

A vector-valued function F ∈ C 1 (Rn , Rm ). is said to be Fréchet differentiable at

x ∈ Rn . if there exists a bounded linear operator A ∈ C(Rn , Rm ). such that

|F (x + h) − F (x) − Ah|
. lim = 0.
h→0 h

We may regard the bounded linear operator A ∈ C(Rn , Rm ). as a matrix in Rm×n .

and write F (x) = AT ∈ Rn×m .. From the definition, we obtain
1.5 Fréchet Derivative 13

⎛ ⎞
∂F1 ∂Fm
∂x1 (x) ··· ∂x1 (x)
⎜ ⎟
F (x) = ⎜
.
⎝
..
.
..
.
⎟ ∈ Rn×m .
⎠ (1.17)
∂F1 ∂Fm
∂xn (x) · · · ∂xn (x)

If m = n., we call A = [F (x)]T . the Jacobian matrix of F at x ∈ Rn ..

If f ∈ C 2 (Rn , R)., we define the Hessian matrix f (x) ∈ Rn×n . as the Fréchet
derivative of f ∈ C 1 (Rn , Rn ).. In matrix form,
⎛ ∂2f ∂2f
⎞
∂x1 ∂x1 (x) ··· ∂x1 ∂xn (x)
⎜ ⎟
f (x) = ⎜
.
⎝
..
.
⎟ ∈ Rn×n .
..
⎠
. (1.18)
∂2f ∂2f
∂xn ∂x1 (x) · · · ∂xn ∂xn (x)

The following identities extend the product rule and chain rule for Fréchet
derivatives.
Proposition 1.5 Let A ∈ Rn×n ., and define f (x) = x T Ax . for x ∈ Rn .. Then

f (x) = (A + AT )x ∈ Rn , f (x) = A + AT .
. (1.19)

Let F, G ∈ C 1 (Rn , Rm )., and define f (x) = F (x)T G(x). for x ∈ Rn .. Then

f (x) = F (x)G(x) + G (x)F (x).

. (1.20)

Given A ∈ Rm×p ., F ∈ C 1 (Rn , Rm )., and G ∈ C 1 (Rn , Rp )., define f (x) =

F (x)T AG(x). for x ∈ Rn .. We have

f (x) = F (x)AG(x) + G (x)AT F (x).

. (1.21)

Given F ∈ C 1 (Rm , Rp ). and G ∈ C 1 (Rn , Rm )., define H (x) = F (G(x)). for

x ∈ Rn .. We have

H (x) = G (x)F (G(x)).

. (1.22)

We are now ready to use Fréchet derivatives to state Taylor expansion in Rn ..

Proposition 1.6 Assume that f ∈ C 2 (Rn , R).. For any u, v ∈ Rn ., there exists w ∈
Rn . such that
1
f (u) = f (v) + (u − v)T f (v) + (u − v)T f (w)(u − v).
. (1.23)
2
Proof Define

g(t) = f (v + (u − v)t), t ∈ [0, 1].

.
14 1 Matrix

By the chain rule in Proposition 1.5,

g (t) = (u − v)T f (v + (u − v)t),

.

and

.g (t) = (u − v)T f (v + (u − v)t)(u − v).

Since f ∈ C 2 (Rn , R)., we have g ∈ C 2 (R, R).. By the one-dimensional Taylor’s

theorem, there exists ξ ∈ (0, 1). such that

1
g(1) = g(0) + g (0) + g (ξ ).
.
2

Setting w = v + (u − v)ξ ∈ Rn . gives (1.23), completing the proof.

Problems
1.1 For any A ∈ Rn×m . and B ∈ Rm×p ., prove that (AB)T = B T AT ..
1.2 If A, B ∈ Rn×n . are invertible, prove that AB is also invertible with (AB)−1 =
B −1 A−1 ..
1.3 Define f (A) = A11 A22 − A21 A12 . for A ∈ R2×2 .. Verify that f is a normalized,
antisymmetric, and multilinear function on R2×2 ..
1.4 Show that a multilinear function f on Rn×n . is antisymmetric if and only if

f (A1 , · · · , An ) = 0
.

for any A1 , · · · , An ∈ Rn . with Aj = Ak . for some j k ..

1.5 Prove Proposition 1.1.
1.6 Given θ ∈ [0, π/2]., evaluate det(A)., where
⎛ ⎞
2 cos θ 1
⎜ 1 2 cos θ 1 ⎟
⎜ ⎟
⎜ .. .. .. ⎟
.A = ⎜ . . . ⎟ ∈ Rn×n .
⎜ ⎟
⎝ 1 2 cos θ 1 ⎠
1 2 cos θ
⎧
⎪
⎪ i = j,
⎨2 cos θ,
Note that A is a tridiagonal matrix with Aij = 1, |i − j | = 1, .
⎪
⎪
⎩0, |i − j | > 1.
1.5 Fréchet Derivative 15

1.7 Evaluate det(B)., where

⎛ ⎞
1 1 ··· 1
⎜ x1 x2 · · · xn ⎟
⎜ ⎟
B=⎜
. .. .. .. ⎟ ∈ R .
n×n
⎝ . . . ⎠
x1n−1 x2n−1 · · · xnn−1

Note that B is a Vandermonde matrix with Bij = xji−1 ..

1.8 For any A ∈ Rm×n . and B ∈ Rn×m ., prove that det(I − AB) = det(I − BA)..
1.9 Prove Proposition 1.5.
Chapter 2
Linear Regression

Abstract In this chapter, we study the linear regression model as a fundamental

example and a cornerstone of machine learning. Two equivalent approaches—
maximum likelihood estimation and least squares approximation—are used to
formulate the associated optimization problem. We derive the optimal solution and
discuss its key properties. In addition, we introduce the variance inflation factor as
a quantitative measure of linear correlation among features.

2.1 Notations

Throughout this chapter, we use d to denote the dimension (i.e., number of features)
of the input data. The data size is n. The output data is a column vector
⎛ ⎞
y1
⎜y2 ⎟
⎜ ⎟
.y = ⎜ . ⎟ ∈ R .
n
(2.1)
⎝ .. ⎠
yn

To account for the constant term in the linear regression model, we introduce an
additional column of ones and denote the augmented input data as
⎛ ⎞
x10 x11 · · · x1d
⎜x20 x21 · · · x2d ⎟
⎜ ⎟
.x = ⎜ . .. ⎟ ∈ R
n×(d+1)
.. , (2.2)
⎝ .. . . ⎠
xn0 xn1 · · · xnd

where xi0 = 1. for i = 1, · · · , n.. For each j = 0, · · · , d ., we denote the j -th column
of x by

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 17

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
18 2 Linear Regression

⎛ ⎞
x1j
⎜x2j ⎟
⎜ ⎟
.xj = ⎜ . ⎟ ∈ R .
n
(2.3)
⎝ .. ⎠
xnj

In particular,
⎛ ⎞
1
⎜1⎟
⎜ ⎟
.x0 = ⎜ . ⎟ ∈ R
n
(2.4)
⎝ .. ⎠
1

is the column vector added to the input data to represent the constant term in the
linear regression model. The other column vectors x1 , · · · , xd ∈ Rn . are also called
the features of the input data, while the output vector y ∈ Rn . is referred to as the
observation.

2.2 Linear Model

The main problem in linear regression model is stated as follows, with illustration
given in Fig. 2.1.

Objective
Given features x1 , · · · , xd ∈ Rn . and observations y ∈ Rn ., we need to find a
parameter vector
⎛ ⎞
p0
⎜ p1 ⎟
⎜ ⎟
.p = ⎜ . ⎟ ∈ R
d+1
(2.5)
⎝ .. ⎠
pd

such that the prediction vector

Y := xp = x0 p0 + x1 p1 + · · · + xd pd
. (2.6)

is a good approximation of the observation vector y, where x =

(x0 , x1 , · · · , xd ) ∈ Rn×(d+1) . is the input data, see (2.2), and x0 ∈ Rn . is
defined as in (2.4).
2.3 Maximum Likelihood Estimation 19

Fig. 2.1 Linear regression model of size n and dimension d

How to determine whether an approximation is good? There are two interpreta-

tions, which are essentially equivalent.

2.3 Maximum Likelihood Estimation

We assume that the error vector

ε := y − Y ∈ Rn
. (2.7)

satisfies a multivariate normal distribution with mean 0 ∈ Rn . and covariance Σ =

σ 2 I ∈ Rn×n ., where I stands for the identity matrix of dimension n. In other words,
if
d
.εi := yi − Yi = yi − (xp)i = yi − xij pj (2.8)
j =0

denotes the i-th component of the error vector (with i = 1, · · · , n.), then ε1 , · · · , εn .
are independent and identically distributed (i.i.d.) random variables satisfying a
normal distribution with mean 0 and variance σ 2 ..
For a given observation vector y ∈ Rn ., the likelihood function is
⎧ ⎛ ⎞2 ⎫
⎪
⎨ 1 n d ⎪
⎬
1 ⎝yi − ⎠
L(p, σ 2 |y) =
. exp − 2 xij pj . (2.9)
(2π σ 2 )n/2 ⎪
⎩ 2σ ⎪
⎭
i=1 j =0

We then define a good approximation as the one that maximizes the above
likelihood function. The problem of finding the maximum likelihood estimation
of the parameter vector p is equivalent to the problem of minimizing the sum of
squared errors (SSE) defined as below:
20 2 Linear Regression

⎛ ⎞2
n d
SSE :=
. ⎝yi − xij pj ⎠ y−Y 2
2, (2.10)
i=1 j =0

where 2 . denotes the l 2 .-norm in Rn ..

2.4 Least Squares Approximation

Another interpretation of the good approximation is the one that minimizes the l 2 .-
norm of the error vector y − Y ..

Optimization Problem
Let x ∈ Rn×(d+1) . and y ∈ Rn . be given. Find p ∈ Rd+1 . such that the objective
function
⎛ ⎞2
n d
1 1 ⎝yi −
.J (p) := y − xp 2
2 = xij pj ⎠ (2.11)
2 2
i=1 j =0

is minimized.

The additional factor 1/2. in the objective function does not affect the solution
of the optimization problem and is included only for the sake of convenience.
The solution of the above optimization problem is also called the least squares
approximation. Since the objective function J (p). in (2.11) and the sum of squared
errors (SSE) in (2.10) differ by a constant factor, the least squares approximation is
equivalent to the maximum likelihood approximation.
Theorem 2.1 If det(x T x) 0., then the objective function J (p). in (2.11) has a
unique global minimum at

p = (x T x)−1 (x T y);
. (2.12)

namely, J (q) > J (p). for any q ∈ Rd+1 . and q p..

Proof Let p be given as in (2.12). For any q ∈ Rd+1 . with q p., we consider the
function

f (t) := J (tq + (1 − t)p),

.
2.5 Sum of Squared Errors 21

for t ∈ R.. Obviously, f (0) = J (p). and f (1) = J (q).. By Taylor expansion, we
obtain
1
f (1) = f (0) + f (0) + f (s),
.
2
for some s ∈ (0, 1).. Note that

f (0) = (q − p)T J (p) = (q − p)T (x T xp − x T y) = 0,

.

and

f (s) = (q −p)T J (sq +(1−s)p)(q −p) = (q −p)T x T x(q −p)

. x(q −p) 22 .

We obtain
1 1
J (q) − J (p) = f (1) − f (0) =
. f (s) = x(q − p) 2
2 ≥ 0.
2 2

Since q − p 0. and det(x T x) 0., we have J (q) > J (p).. Otherwise, x(q −
p) = 0. and x T x(q−p) = 0.; namely, 0 is an eigenvalue of x T x . with a corresponding
eigenvector q − p ., which is a contradiction to the assumption det(x T x) 0.. This
completes the proof.

2.5 Sum of Squared Errors

In this section, we assume det(x T x) 0., and let p be the optimal solution as given
in (2.12). The following proposition shows that the prediction vector Y = xp. is the
projection of the observation vector y on the linear subspace spanned by the vectors
x0 , · · · , xd .; see Fig. 2.2.
Proposition 2.1 Given x ∈ Rn×(d+1) . and y ∈ Rn . with det(x T x) 0., let Y = xp .
with p = (x T x)−1 (x T y).. We have the following extension of Pythagorean theorem

. y−Y 2
2 y 2
2 Y 2
2. (2.13)

Proof Note that

Y T (y − Y ) = pT x T y − p T x T xp = pT (x T y − x T xp) = 0.
.

We have Y T Y = Y T y = y T Y ., and hence,

. y−Y 2
2 = yT y − Y T y − yT Y + Y T Y = yT y − Y T Y y 2
2 Y 2
2.
22 2 Linear Regression

Fig. 2.2 Pythagorean theorem

This completes the proof.

Denote

xT x xT y
B := (x, y)T (x, y) =
. ∈ R(d+2)×(d+2) . (2.14)
yT x yT y

The following proposition shows that the last diagonal term (i.e., the (d +2, d +2).
entry) of B −1 . is the reciprocal of the SSE.
Proposition 2.2 Given x ∈ Rn×(d+1) . and y ∈ Rn . with det(x T x) 0., let Y = xp .
with p = (x T x)−1 (x T y).. We have

−1
xT x xT y (x T x)−1 + pλpT − pλ
. = , (2.15)
yT x yT y −λpT λ

where
1 1
λ=
. = .
SSE y−Y 2
2

Proof It is easily seen from x T xp = x T y . that

xT x xT y xT x 0 I p
. = .
yT x yT y 0 1 yT x yT y

Assume

I p A α I 0
. = .
yT x yT y βT λ 0 1
2.6 Variance Inflation Factor 23

We then have

A + pβ T = I, y T xA + y T yβ T = 0,
.

α + pλ = 0, y T xα + y T yλ = 1.

Solving these four equations yields

1 1 1
α = −pλ, λ =
. = T = ,
yT y − y xp
T y (y − Y ) SSE

and

.β T = −λy T x, A = I + pλy T x.

Finally,

−1
xT x xT y I + pλy T x − pλ (x T x)−1 0
. =
yT x yT y −λy T x λ 0 1

(x T x)−1 + pλy T x(x T x)−1 − pλ

=
−λy T x(x T x)−1 λ

(x T x)−1 + pλpT − pλ
= .
−λpT λ

This completes the proof.

2.6 Variance Inflation Factor

The sample variance of the output data is defined as

n
1
σy2 :=
. (yi − ȳ)2 , (2.16)
n
i=1

where
n
1
ȳ :=
. yi (2.17)
n
i=1

is the mean/average of the output data. Assume the input data x ∈ Rn×(d+1) . satisfies
the condition det(x T x) 0.. Assume further that the sample variance of the output
24 2 Linear Regression

data is positive. We then define the coefficient of determination (R squared) of the

linear regression model as
n
i=1 (yi − Yi )
SSE 2
. R 2 := 1 − =1− n , (2.18)
nσy2 i=1 (yi − ȳ)
2

where SSE y−Y 2 . and

2 Y = xp. with p = (x T x)−1 (x T y)..
Theorem 2.2 The coefficient of determination is bounded by 1 and nonnegative:
0 ≤ R 2 ≤ 1.. Moreover, R 2 = 1. if and only if the SSE is zero.
Proof We only need to show that R 2 ≥ 0.. The other conclusions are obvious.
Denote q = (ȳ, 0, · · · , 0)T ∈ Rd+1 .. Recall the definition of the objective function
in (2.11). It then follows from (2.18) that

J (p)
R2 = 1 −
. .
J (q)

By Theorem 2.1, we obtain J (p) ≤ J (q)., and hence R 2 ≥ 0.. This completes the
proof.
Recall that x = (x0 , x1 , · · · , xd ) ∈ Rn×(d+1) . with x0 = (1, 1, · · · , 1)T ∈
Rn .. If R 2 = 1., then the observation vector y can be expressed as a linear
combination of x0 , · · · , xd ., and the matrix B = (x, y)T (x, y). becomes singular.
For each k = 1, · · · , d ., we consider the linear regression model with x (k) :=
(x0 , · · · , xk−1 , xk+1 , · · · , xd ) ∈ Rn×d . as the input data and xk . as the output data;
that is, we find the parameter vector p (k) ∈ Rd . such that the prediction vector
Xk = x (k) p(k) ∈ Rn . has the smallest Euclidean distance from xk ∈ Rn .. The
corresponding sum of squared errors is denoted as
n
SSEk
. xk − Xk 2
2 = (xik − Xik )2 . (2.19)
i=1

The sample variance of xk . is given by

n
1
2
.σk := (xik − x̄k )2 , (2.20)
n
i=1

where
n
1
x̄k :=
. xik (2.21)
n
i=1

is the mean/average of xk .. Similar to the original linear regression model, we define

the coefficient of determination (R squared) as
2.7 Python Code 25

n
i=1 (xik − Xik )
SSEk 2
Rk2 := 1 −
. =1− n . (2.22)
nσk2 i=1 (xik − x̄k )
2

Finally, the variance inflation factor (VIF) is defined as

n
nσk2 i=1 (xik − x̄k )
1 2
V I Fk =
. = = n . (2.23)
1 − Rk2 i=1 (xik − Xik )
SSEk 2

The following theorem indicates that a large VIF implies a strong linear
correlation between the features.
Theorem 2.3 Given the input data x = (x0 , x1 , · · · , xd ) ∈ Rn×(d+1) . with x0 =
(1, 1, · · · , 1)T ∈ Rn ., let V I F1 , · · · , V I Fd . be the variance inflation factors defined
as in (2.23). We obtain det(x T x) = 0. if and only if V I Fk = ∞. for some k ∈
{1, · · · , d}..
Proof Note that det(x T x) = 0. if and only if there exists k ∈ {1, · · · , d}. such that
xk . is linearly dependent of x0 , · · · , xk−1 , xk+1 , · · · , xd ., which is also equivalent to
SSEk = 0., where SSEk . is the sum of squared errors defined in (2.19). On account
of (2.23), we have SSEk = 0. if and only if V I Fk = ∞.. This completes the proof.

The following proposition provides an alternative formula of the VIF.

Proposition 2.3 Given the input data x = (x0 , x1 , · · · , xd ) ∈ Rn×(d+1) . with x0 =
(1, 1, · · · , 1)T ∈ Rn ., assume det(x T x) 0., and let λ0 , · · · , λd . be the diagonal
terms of (x T x)−1 .. For each k = 1, · · · , d ., let V I Fk . be the variance inflation factor
defined as in (2.23). We have V I Fk = nσk2 λk ., where σk2 . is the sample variance of
xk . as given in (2.20).
Proof By permutation, we may assume without loss of generality that k = d .. A
simple application of Proposition 2.2 gives λk = 1/SSEk ., where SSEk . is the sum of
squared errors defined in (2.19). In view of (2.23), we obtain V I Fk = nσk2 /SSEk =
nσk2 λk .. This completes the proof.

2.7 Python Code

In this section, we provide a Python code to find the least squares solution,
coefficient of determination, and variance inflation factors of the linear regression
model.
26 2 Linear Regression

Linear Regression Code

import numpy as np

def linear_regression(x,y):
n=[Link][0]
p=[Link]([Link]([Link](x.T,x)),[Link](x.T,y))
SSE=[Link]([Link](x,p))**2
R2=1-SSE/(n*[Link](y))
Lambda=[Link]([Link]([Link](x.T,x)))
VIF=[Link](Lambda,n*[Link](x,axis=0))[1:]
result={
’p’:p,
’R2’:R2,
’VIF’:VIF
}
return result

n=1000
d=4
x0=[Link]((n,1))
x=[Link]((x0,[Link]((n,d))),axis=1)
p_exact=[Link]((d+1,1))
y=[Link](x,p_exact)+[Link](0,0.1,(n,1))

result=linear_regression(x,y)
p=result[’p’]
R2=result[’R2’]
VIF=result[’VIF’]
print(R2)
print(VIF)
print([Link]((p,p_exact),axis=1))

In the code, we choose the sample size n = 1000. and the dimension d = 4.. The
input data x1 , · · · , xd ∈ Rn . is randomly generated from the uniform distribution.
The exact parameter vector is set to be (1, · · · , 1)T ∈ Rd+1 .. The output data is
perturbed by normal random variables. A sample result is given below:
2.7 Python Code 27

Sample Result
0.9702959227690222
[1.00333685 1.005423 1.00388619 1.00681191]
[[0.99675087 1. ]
[1.00069386 1. ]
[0.98648815 1. ]
[1.00547332 1. ]
[1.01346316 1. ]]

It is noted that the coefficient of determination is very close to 1, the variance

inflation factors are very close to 1, and the least squares solution is also very close
to the exact solution.
To compute p = (x T x)−1 (x T y)., we adopt a code that, while being intuitively
straightforward, regrettably suffers from low efficiency. Nonetheless, for data with
a comparatively small dimension d, the utilization of this code is still tolerable. It
must be stressed, however, that in actual practical applications, whenever feasible,
one should avoid directly computing the inverse of x T x .. This is because such inverse
computations generally involve substantial computational costs and complexity,
inevitably resulting in poor efficiency. Hence, during real-world operations, it is
essential to seek out more optimized computational pathways. In Chap. 9, we will
introduce a collection of gradient methods to approximate the solution for the
linear system x T xp = x T y ., providing an alternative to circumvent the potential
inefficiencies linked with direct computation of (x T x)−1 . with large dimension d.

Problems
2.1 For any x ∈ Rn×(d+1) ., show that x T x . is a symmetric and positive semidefinite
matrix in R(d+1)×(d+1) .. If x T x . is invertible, then it is positive definite.
2.2 Let x = (x0 , x1 , · · · , xd ) ∈ Rn×(d+1) ., where x0 , · · · , xd ∈ Rn . and x0 0..
Prove that det(x T x) = 0. if and only if there exists k ∈ {1, · · · , d}. such that xk . can
be represented as a linear combination of x0 , · · · , xk−1 , xk+1 , · · · , xd .; namely,

xk = c0 x0 + · · · + ck−1 xk−1 + ck+1 xk+1 + · · · + cd xd ,

.

for some c0 , · · · , ck−1 , ck+1 , · · · , cd ∈ R..

2.3 Given x ∈ Rn×(d+1) . and y ∈ Rn ., define a function
⎛ ⎞2
n d
1 1 ⎝yi −
.J (p) = y − xp 2
2 = xij pj ⎠
2 2
i=1 j =0
28 2 Linear Regression

for
⎛ ⎞
p0
⎜ p1 ⎟
⎜ ⎟
.p = ⎜ . ⎟ ∈ R
d+1
.
⎝ .. ⎠
pd

(a) Verify that J (p) = x T (xp − y) ∈ Rd+1 .; namely,

⎛ ⎞
n d
∂
. J (p) = [J (p)]k = [x T (xp − y)]k = xik ⎝ xij pj − yi ⎠ ,
∂pk
i=1 j =0

for k = 0, · · · , d ..
(b) Verify that J (p) = x T x ∈ R(d+1)×(d+1) .; namely,
n
∂2
. J (p) = [J (p)]j k = (x T x)j k = xij xik ,
∂pj ∂pk
i=1

for j, k = 0, · · · , d ..
2.4 Given J ∈ C 2 (Rd+1 , R)., for any p ∈ Rd+1 . and q ∈ Rd+1 ., prove that there
exists ξ ∈ Rd+1 . such that

1
.J (q) = J (p) + (q − p)T J (p) + (q − p)T J (ξ )(q − p).
2
2.5 Prove that the variance inflation factor is no less than one; namely, V I Fk ≥ 1.,
where V I Fk . is defined in (2.23).
2.6 Given the input data x = (x0 , x1 , · · · , xd ) ∈ Rn×(d+1) . with x0 =
(1, 1, · · · , 1)T ∈ Rn ., let SSEk . with k = 1, · · · , d . be the sum of squared errors
defined as in (2.23). Assume det(x T x) 0., and let λ0 , · · · , λd . be the diagonal
terms of (x T x)−1 .. Show that λk = 1/SSEk . for k = 1, · · · , d ..
2.7 Use a linear regression model to solve the problem of fitting the points
(0, 0), (1, 2), (2, 1), (3, 3). by a straight line.
2.8 Given a positive integer n, calculate the mean and variance of the sequence
1, 2, · · · , n..

2.9 Given an input matrix

⎛ ⎞
110
.x = ⎝1 0 1⎠ ,

100
2.7 Python Code 29

calculate the variance inflation factors as defined in (2.23). Note that the sample size
is n = 3. and the dimension is d = 2..
2.10 Given an input matrix x ∈ Rn×(d+1) . such that det(x T x) 0., fix an
observation vector y ∈ Rn ., and find the parameter vector p ∈ Rd+1 . and the variance
σ 2 > 0. such that the likelihood function L(p, σ 2 |y). defined in (2.9) is maximized.
Chapter 3
Regularization

Abstract A large variance inflation factor in a linear regression model indicates

strong linear correlation among features and may lead to overfitting. In this chapter,
we address this issue by regularizing the objective function with the l 2 .-norm of
the parameter vector and introduce the ridge regression model. We also study the
LASSO (least absolute shrinkage and selection operator) method, which employs
l 1 .-regularization. Since the l 1 .-norm is not differentiable, we introduce the concept
of the subdifferential and formulate the problem within the framework of convex
optimization. An iterative algorithm is then derived to compute a numerical solution
for linear regression with l 1 .-regularization.

3.1 Overfitting Problem

We continue with the linear regression code in Sect. 2.7. Now, we assume x2 . and x1 .
are closely related in the sense that the difference x2 − x1 . is sampled from a normal
distribution with mean 0 and a small variance. For simplicity, we do not repeat the
code for the linear regression function in Sect. 2.7.

Overfitting Example
import numpy as np

n=1000
d=4
x0=[Link]((n,1))
x=[Link]((x0,[Link]((n,d))),axis=1)
x[:,2]=x[:,1]+[Link](0,0.01,(n,))
p_exact=[Link]((d+1,1))
y=[Link](x,p_exact)+[Link](0,0.1,(n,1))

(continued)

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 31

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
32 3 Regularization

result=linear_regression(x,y)
p=result[’p’]
R2=result[’R2’]
VIF=result[’VIF’]
print(R2)
print(VIF)
print([Link]((p,p_exact),axis=1))

A sample result is given below:

Overfitting Result
0.9796829459418794
[803.11874709 802.97720154 1.01312437 1.01804729]
[[1.00094331 1. ]
[1.19645218 1. ]
[0.7916922 1. ]
[1.00670459 1. ]
[1.0056871 1. ]]

It is noted that the coefficient of determination is still close to 1. However, the

variance inflation factors V I F1 . and V I F2 . are large. This indicates that the features
may be linearly correlated and the model may encounter a problem of overfitting.
The least squares approximation shows that p1 . and p2 . are not close to 1, while the
sum p1 + p2 . is very close to 1 + 1 = 2..
In practice, when we deal with a high-dimension problem with many features,
it is very likely that some features are linearly correlated. In this case, we may find
many solutions with small bias (i.e., the SSE is small). In particular, if det(x T x) =
0., then the linear regression model has infinitely many solutions. Moreover, when
the variance is high (i.e., the VIF is large), the condition number of x T x . becomes
large, and the solution of x T xp = x T y . has low accuracy. To balance the bias
and variance, we need to introduce the l 2 .-regularization and consider the ridge
regression problem.
3.2 Ridge Regression 33

3.2 Ridge Regression

In the case of large VIF (i.e., x T x . is ill-conditioned), we need to modify the

optimization problem in Sect. 2.4 with the following objective function:

1 λ
J (p) =
. y − xp 2
2 + p 22 , (3.1)
2 2
where λ > 0. is a regularization constant to be tuned. It is also called a
hyperparameter or a tuning parameter.
Theorem 3.1 For any x ∈ Rn×(d+1) ., y ∈ Rn ., and λ > 0., the objective function
J (p). in (3.1) has a unique global minimum at

p = (x T x + λI )−1 (x T y);
. (3.2)

namely, J (q) > J (p). for any q ∈ Rd+1 . and q p..

Proof Let p be given as in (3.2). We have

J (p) = x T (xp − y) + pλ = (x T x + λI )p − x T y = 0.
.

For any q ∈ Rd+1 . with q p., we obtain from Taylor expansion

1 λ
J (q) − J (p) =
. (q − p)T (x T x + λI )(q − p) ≥ (q − p) 2
2 > 0.
2 2
This completes the proof.
We come back to the overfitting problem in Sect. 3.1. The Python code imple-
menting the ridge regression is given below:

Ridge Regression
import numpy as np

n=1000
d=4
x0=[Link]((n,1))
x=[Link]((x0,[Link]((n,d))),axis=1)
x[:,2]=x[:,1]+[Link](0,0.01,(n,))
p_exact=[Link]((d+1,1))
y=[Link](x,p_exact)+[Link](0,0.1,(n,1))

(continued)
34 3 Regularization

p=[Link]([Link]([Link](x.T,x)),[Link](x.T,y))
la=1
A=la*[Link](d+1)+[Link](x.T,x)
p_ridge=[Link]([Link](A),[Link](x.T,y))
print([Link]((p,p_exact,p_ridge),axis=1))

We choose the regularization parameter λ = 1.. A sample result is given below:

Sample Result
[[0.98693898 1. 0.99435015]
[1.15005244 1. 1.01000321]
[0.86133075 1. 0.99605137]
[1.00696983 1. 1.00073844]
[0.99278865 1. 0.98752897]]

It is noted that the result from the ridge regression (linear regression with l 2 .-
regularization) is closer to the exact parameter vector than that of linear regression.
This example illustrates that l 2 .-regularization helps reduce the side effect of
overfitting the data with correlated features.

3.3 Convex Optimization

Another commonly used regularization method is adding the l 1 .-norm of the

parameter vector to the objective function. As we shall see later, this regularization
helps reduce the number of nonzero parameters in the estimation and enables us to
select the significant features of the input data. The method is also called LASSO,
which stands for least absolute shrinkage and selection operator.
Since l 1 .-norm may not be differentiable at some points, we need to generalize
the definitions of derivatives and the ideas of critical points. First, we recall the
definition of a convex function in an interval. Let I ⊂ R. be an open interval. A
continuous function f ∈ C(I, R). is convex if and only if

f (tu + (1 − t)v) ≤ tf (u) + (1 − t)f (v),

. (3.3)

for all u ∈ I ., v ∈ I ., and t ∈ [0, 1].. Note that a convex function is not necessarily
differentiable. For example, f (u) = |u|. is convex in R. but not differentiable at
u = 0.. The following lemma gives an alternative definition of a convex function;
see Fig. 3.1.
3.3 Convex Optimization 35

Fig. 3.1 Convex function

Lemma 3.1 Let I ⊂ R. be an open interval and f ∈ C(I, R).. Then f is convex in
I if and only if for any u ∈ I ., v ∈ I ., and w ∈ I . such that u < w < v ., we have

f (w) − f (u) f (v) − f (u) f (v) − f (w)

. ≤ ≤ . (3.4)
w−u v−u v−w

Proof First, we assume that (3.4) holds for all u ∈ I ., v ∈ I ., and w ∈ I . such that
u < w < v .. In particular, by choosing w = tu + (1 − t)v . with t ∈ (0, 1)., we obtain
from the first inequality in (3.4) that

f (w) − f (u) f (v) − f (u)

. ≤ .
(1 − t)(v − u) v−u

A further simplification yields (3.3). When t = 0. or t = 1., the inequality (3.3)

becomes an obvious equality.
Next, we assume that (3.3) holds for all u ∈ I ., v ∈ I ., and t ∈ [0, 1].. Given any
w ∈ (u, v)., we choose t = (v −w)/(v −u) ∈ (0, 1). such that w = tu+(1−t)v .. The
inequalities in (3.4) can be easily obtained from (3.3). This completes the proof.
Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R.. For any u ∈ I .,
one can prove that the one-sided derivatives

f (u + h) − f (u)
f ± (u) := lim
.
h→0± h
36 3 Regularization

exist and f − (u) ≤ f + (u).. We then define the subdifferential of f at u as the

following closed interval:

∂f (u) := [f − (u), f + (u)].

. (3.5)

Any value in ∂f (u). is called a subderivative of f at u. If f is differentiable

at u, then the subdifferential shrinks to a single point f (u)., and the subderivative
coincides with the derivative. It is well known that if f is convex and differentiable
on an open interval I ⊂ R., then u ∈ I . is a local minimum of f if and only if
f (u) = 0.. Moreover, any local minimum of f in I is global in I . By using the
subdifferential, we can extend this result in the following proposition.
Proposition 3.1 Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R..
We have the following results:
(a) u ∈ I . is a local minimum point of f if and only if 0 ∈ ∂f (u)..
(b) If u ∈ I . and v ∈ I . with u < v . are two local minimum points of f , then
f (w) = 0. and f (w) = f (u) = f (v). for any w ∈ (u, v)..
Proof If u ∈ I . is a local minimum point of f , then f (u + h) ≥ f (u). for any
sufficiently small h > 0., and hence f + (u) ≥ 0.. Similarly, we have f − (u) ≤ 0..
Therefore, 0 ∈ ∂f (u).. On the other hand, if 0 ∈ ∂f (u)., namely, f − (u) ≤ 0 ≤
f + (u)., then we have

. lim g(h) ≥ 0,
h→0+

where g(h) := [f (u + h) − f (u)]/ h. is the difference quotient. By Lemma 3.1,

g(h1 ) ≥ g(h2 ). for any small h1 > h2 > 0.. Hence, g(h) ≥ 0. and f (u + h) ≥ f (u).
for all small h > 0.. Similarly, we can show that f (u − h) ≥ f (u). for all small
h > 0.. Therefore, u is a local minimum of f . This proves (a).
Now, we assume that u ∈ I . and v ∈ I . with u < v . are two local minimum points
of f . By (a), we have

f + (u) ≥ 0 ≥ f − (v).
.

By letting w → u+ . in the first inequality of (3.4) and w → v − . in the second

inequality of (3.4), respectively, we obtain

f (v) − f (u)
f + (u) ≤
. ≤ f − (v).
v−u

Coupling the above two inequalities gives

.f + (u) = f − (v) = f (v) − f (u) = 0.

3.3 Convex Optimization 37

Moreover, for any w ∈ (u, v)., we obtain from

0 = f + (u) ≤ f − (w) ≤ f + (w) ≤ f − (v) = 0

.

that f (w) = 0.. A simple integration from u to w then gives f (w) = f (u) = f (v).
for all w ∈ (u, v).. This proves (b).
Denote f (u) := λ|u|.. A simple calculation gives
⎧
⎪
⎪ u > 0,
⎨λ,
∂f (u) =
. −λ, u < 0, (3.6)
⎪
⎪
⎩[−λ, λ], u = 0.

Consider a simple optimization problem with the objective function

1
J (u) =
. (u − z)2 + λ|u|, u ∈ R, (3.7)
2
where z ∈ R. and λ > 0. are given. Denote
⎧
⎪
⎨z − λ,
⎪ z > λ,
Fλ (z) :=
. z + λ, z < −λ, (3.8)
⎪
⎪
⎩0, z ∈ [−λ, λ].

Lemma 3.2 The objective function J (u). in (3.7) is minimized at u∗ = Fλ (z). as

defined in (3.8).
Proof By Proposition 3.1, u∗ ∈ R. is a local minimum of J (u). if and only if z−u∗ ∈
∂f (u∗ ).. In view of (3.6), we obtain (i) z − u∗ = λ. if u∗ > 0., (ii) z − u∗ = −λ. if
u∗ < 0., and (iii) z − u∗ ∈ [−λ, λ]. if u∗ = 0.. This proves u∗ = Fλ (z)..
The following theorem shows that the local minimum of a sum of f (u) = λ|u|.
with λ > 0. and a differentiable function g can be interpreted as a fixed point of a
function involving Fλ . defined in (3.8).
Theorem 3.2 Let g ∈ C 1 (I, R). be differentiable in an open interval I ⊂ R.. Set
f (u) = λ|u|. with λ > 0.. Then u∗ ∈ I . is a local minimum of f (u) + g(u). if and
only if u∗ = Fλα (u∗ − αg (u∗ )). for any α > 0.. Here, Fλ . is defined in (3.8).
Proof Note that u∗ ∈ I . is a local minimum of f (u) + g(u). if and only if 0 ∈
∂f (u∗ ) + g (u∗ ).. For any α > 0., we also note that u∗ ∈ I . is a local minimum
of αf (u) + [u − u∗ + αg (u∗ )]2 /2. if and only if 0 ∈ α∂f (u∗ ) + αg (u∗ ).. Hence,
u∗ ∈ I . is a local minimum of f (u) + g(u). if and only if it is a local minimum
of αf (u) + [u − u∗ + αg (u∗ )]2 /2. for any α > 0., which satisfies the fixed point
equation u∗ = Fλα (u∗ − αg (u∗ )).. This completes the proof.
38 3 Regularization

3.4 LASSO

We are now ready to investigate the optimization problem with l 1 .-regularization,

which is also referred to as LASSO (least absolute shrinkage and selection
operator).

LASSO
Let x ∈ Rn×(d+1) . and y ∈ Rn . be given. Find p ∈ Rd+1 . such that the objective
function
⎛ ⎞2
n d d
1 1 ⎝yi −
J (p) :=
. y − xp 2
2 +λ p 1 = xij pj ⎠ + λ |pj |
2 2
i=1 j =0 j =0
(3.9)
is minimized.

We shall extend Theorem 3.2 from one-dimensional real line to a higher

dimensional space.
Theorem 3.3 Let Fλ . be defined as in (3.8). Given x ∈ Rn×(d+1) . and y ∈ Rn ., fix
any α > 0.. Then p ∈ Rd+1 . is a local minimum of the objective function J (p). in
(3.9) if and only if

pj = Fλα (pj − α(x T xp − x T y)j ),

.

for j = 0, · · · , d ..
Proof Define

1
g(p) :=
. y − xp 22 .
2

Note that g (p) = x T (xp − y).. For each j = 0, · · · , d ., we have

n d n
∂
. g(p) = (x T xp − x T y)j = xij xik pk − xij yi .
∂pj
i=1 k=0 i=1

By regarding g(p). as a one-variable function of pj ., we obtain from Theorem 3.2

our desired formula. This completes the proof.
For simplicity, we extend the definition of Fλ . in (3.8) to be a function in
C(Rd+1 , Rd+1 ). such that
3.5 Python Code 39

[Fλ (u)]j = Fλ (uj )

. (3.10)

for u = (u0 , u1 , · · · , ud ) ∈ Rd+1 . and j = 0, 1, · · · , d .. Theorem 3.3 provides an

iterative method of finding a local minimum of J (p). in (3.9).

p(k+1) = Fλα (p(k) − αx T (xp(k) − y)).

. (3.11)

Here, λ > 0. is the regularization parameter, and α > 0. is also called the learning
rate. Both hyperparameters are to be tuned in practical applications. In the special
case λ → 0+ ., the above iteration reduces to the steepest descent method:

p(k+1) = p(k) − αx T (xp(k) − y) = p(k) − αJ (p(k) ).

.

Comparing with the direct formula in (2.12), the steepest descent method avoids
solving the inverse of the matrix x T x . which is generally dense and sometimes ill-
conditioned. The choice of the learning rate α > 0. depends on the specific problem.
In principle, a large α . may lead to divergence of the iteration, while a small α . may
induce some computation costs.

3.5 Python Code

In this section, we provide a Python code to implement the iterative algorithm (3.11)
in LASSO.

LASSO Code
import numpy as np

def F(u,la):
return (u>la)*(u-la)+(u<-la)*(u+la)

def lasso(x,y,la,al,iter_max):
d=[Link][1]-1
p=[Link]((d+1,1))
for iter in range(iter_max):
u=p-al*[Link](x.T,[Link](x,p)-y)
p=F(u,la*al)
return p

n=1000
(continued)
40 3 Regularization

d=4
x0=[Link]((n,1))
x=[Link]((x0,[Link]((n,d))),axis=1)
p_exact=[Link]((d+1,1))
p_exact[3]=0
p_exact[4]=0
y=[Link](x,p_exact)+[Link](0,0.1,(n,1))

p=[Link]([Link]([Link](x.T,x)),[Link](x.T,y))
p_lasso=lasso(x,y,la=1,al=0.0001,iter_max=10000)
print([Link]((p,p_exact,p_lasso),axis=1))

In the code, we choose the sample size n = 1000. and the dimension d = 4.. The
input data x1 , · · · , xd ∈ Rn . is randomly generated from the uniform distribution.
The exact parameter vector is set to be pexact = (1, 1, 1, 0, 0)T ∈ R5 .. The output
data is perturbed by normal random variables. A sample result is given below:

Sample Result
[[ 1.01146268 1. 1.00860743]
[ 1.00235203 1. 0.99627135]
[ 0.99401163 1. 0.9872306 ]
[-0.0103034 0. 0. ]
[-0.00607228 0. 0. ]]

It is noted that the linear regression (without regularization) gives small but
nonzero estimations to p3 . and p4 .. However, by choosing a suitable regularization
parameter, the LASSO method estimates p3 = 0. and p4 = 0.. This example shows
that the LASSO method is able to select the significant features of the input data.

3.6 Discussions

We use an illustrative example to compare l 2 .-regularization and l 1 .-regularization.

Consider an objective function

1 1
J (p) =
. (p0 − 1)2 + (p1 − 0.01)2 ,
2 2
3.6 Discussions 41

for p = (p0 , p1 )T ∈ R2 .. Clearly, J (p). is minimized at

1
p∗ =
. .
0.01

Now, we fix λ > 0. and consider two objective functions

1 1 λ
Jˆ(p) = (p0 − 1)2 + (p1 − 0.01)2 + (p02 + p12 ),
.
2 2 2
and
1 1
J˜(p) = (p0 − 1)2 + (p1 − 0.01)2 + λ(|p0 | + |p1 |),
.
2 2

respectively, as the l 2 .-regularization and l 1 .-regularization of J (p).. It is easily seen

that Jˆ(p). is minimized at

1 1
p̂ ∗ =
. ,
1 + λ 0.01

which approaches p ∗ . as λ → 0. and (0, 0)T . as λ → ∞.. Finally, we obtain from

Lemma 3.2 that J˜(p). is minimized at

(1 − λ)+
p̃ ∗ =
. , (3.12)
(0.01 − λ)+

where the plus subscript denotes the maximum of the number and zero; namely,

x, x ≥ 0,
x+ := max{x, 0} =
. (3.13)
0, x ≤ 0.

Similar to l 2 .-regularization, we have p̃ ∗ → p∗ . as λ → 0. and p̃ ∗ → (0, 0)T .

as λ → ∞.. An important and useful property of l 1 .-regularization is that only
the significant features are selected if the regularization parameter is appropriately
chosen. In our illustrative example, we have p̃1 = 0. for λ ∈ [0.01, 1).. It is necessary
to point out that p̃ ∗ = (0, 0)T . if λ ≥ 1.. This suggests us that in practice, we should
avoid setting λ. to be too large; otherwise the parameter vector will become zero,
and none of the features will be selected.

Problems
3.1 For any x ∈ Rn×(d+1) . and y ∈ Rn ., show that the equation x T xp = x T y . has
infinitely many solutions if and only if det(x T x) = 0..
42 3 Regularization

3.2 Given the input data

⎛ ⎞
1 0 1
⎜1 1 2⎟⎟
.x = ⎜
⎝1 2 3⎠
1 3 4

and the output data y = (0, 2, 1, 3)T ., find the optimal parameter vector(s) p ∈ R3 .
such that the sum of squared errors SSE y − xp 22 . is minimized.
3.3 Given x ∈ Rn×(d+1) . and λ > 0., prove that the matrix x T x + λI . is positive
definite and hence invertible, where I is the identity matrix of dimension d + 1..
3.4 Given x ∈ Rn×(d+1) ., y ∈ Rn ., and λ > 0., define a function
⎛ ⎞2
n d d
1 λ 1 ⎝yi − 1
.J (p) = y − xp 2+
2
p 2
2 = xij pj ⎠ + pj2
2 2 2 2
i=1 j =0 j =0

for
⎛ ⎞
p0
⎜ ∂1 ⎟
⎜ ⎟
.p = ⎜ . ⎟ ∈ R
d+1
.
⎝ .. ⎠
∂d

(a) Verify that J (p) = x T (xp − y) + pλ ∈ Rd+1 .; namely,

⎛ ⎞
n d
∂
. J (p) = [J (p)]k = xik ⎝ xij pj − yi ⎠ + λpk ,
∂pk
i=1 j =0

for k = 0, · · · , d .. (b) Verify that J (p) = x T x + λI ∈ R(d+1)×(d+1) .; namely,

n
∂2
. J (p) = [J (p)]j k = xij xik + λδj k ,
∂pj ∂pk
i=1

for j, k = 0, · · · , d ., where

1 j = k,
δj k =
.
0 j k

is the Kronecker delta symbol.

3.6 Discussions 43

3.5 In the Python code given in Section 3.2, tune the hyperparameter λ. from 10−3 .
to 103 ., and find the pattern of changes in the bias and variance of the solution.
3.6 Let f ∈ C 2 (I, R)., where I ⊂ R. is an open interval. Prove that f is convex in
I if and only if f (u) ≥ 0. for all u ∈ I ..
3.7 Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R.. For any
u ∈ I ., prove that the one-sided derivatives

f (u + h) − f (u)
f ± (u) := lim
.
h→0± h

exist and f − (u) ≤ f + (u).. Moreover, the derivative f (u). exists if and only if
f − (u) = f + (u)..
3.8 Let f ∈ C 1 (I, R). be a convex function on an open interval I ⊂ R.. Then u ∈ I .
is a local minimum of f if and only if f (u) = 0.. Moreover, any local minimum of
f in I is global in I .
3.9 Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R.. For any
u ∈ I . and v ∈ I . with u < v ., prove that f + (u) ≤ f − (v)..
3.10 Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R.. Given
u ∈ I ., prove that c ∈ ∂f (u). if and only if f (v) − f (u) ≥ c(v − u). for any v ∈ I ..
3.11 Let f ∈ C(I, R). be a convex function on an open interval I ⊂ R.. Given
g ∈ C 1 (I, R)., prove that u ∈ I . is a local minimum of f + g . if and only if 0 ∈
∂f (u) + g (u). (i.e., f − (u) ≤ −g (u) ≤ f + (u).).
3.12 Verify (3.12).
Chapter 4
Nonlinear Regression

Abstract Linear regression is applicable only when the output depends approxi-
mately linearly on the input. In many real-world datasets, however, this relationship
is not necessarily linear, and a nonlinear regression model may provide a more
appropriate description. In this chapter, we use logic gates as illustrative examples to
demonstrate how nonlinear regression succeeds in situations where linear regression
fails. The nonlinear regression model differs from the linear regression model in two
major aspects: (i) The linear activation function is replaced by a sigmoid function (or
a hyperbolic tangent function), and (ii) a logarithmic loss function is used to measure
approximation errors instead of a quadratic loss. The resulting optimization problem
is not well-posed in the sense that a global minimum cannot be attained at any finite
parameter vector. Nevertheless, the gradient descent method can be employed as
an efficient numerical algorithm to obtain a good approximate solution. The logic
gate examples illustrate how machine learning can achieve accurate and efficient
predictions without explicitly solving the open problem of global optimization.

4.1 Nonlinear Data

We use a simple example to illustrate how linear regression fails in fitting nonlinear
data. Consider a logic gate with two binary inputs x1 . and x2 . and one binary output
y. The dimension is d = 2., and the sample size is n = 4.. Table 4.1 illustrates the
input and output data of the AND gate.
For convenience, we use 1 and − 1. to denote True and False, respectively. The
input matrix of the AND gate can be written as
⎛ ⎞
1 1 1
⎜1 1 −1⎟
.x = ⎜ ⎟
⎝1 −1 1 ⎠ ∈ R .
4×3
(4.1)
1 −1 −1

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 45

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
46 4 Nonlinear Regression

Table 4.1 The AND gate x1 . x2 . y

True True True
True False False
False True False
False False False

The output vector of the AND gate is

⎛ ⎞
1
⎜−1⎟
.y = ⎜ ⎟
⎝−1⎠ ∈ R .
4
(4.2)
−1

A simple calculation gives the parameter vector

⎞⎛
−0.5
(x y) = ⎝ 0.5 ⎠ ∈ R3
T −1 T
.p = (x x)

0.5

and the prediction vector

⎛ ⎞
0.5
⎜−0.5⎟
.Y = xp = ⎜ ⎟
⎝−0.5⎠ ∈ R .
4

−1.5

Since the linear regression model fails to give a satisfactory prediction, we need to
use a nonlinear regression model to fit the data.

4.2 Sigmoid Function

The sigmoid function is usually referred to the following function:

1
S(t) :=
. , (4.3)
1 + e−t

which is increasing and bounded by its two horizontal asymptotes 0 and 1. However,
from mathematical point of view, any function obtained by scaling and shifting of
S(t). will have equivalent effects as S(t). in nonlinear regression and hence should
still be called the sigmoid function. In other words, for any real constants α1 ., β1 .,
α0 ., and β0 . such that α1 β1 = 0., the function
4.2 Sigmoid Function 47

β1
u(t) := β1 S(α1 t + α0 ) + β0 =
. + β0 (4.4)
1 + e−α1 t−α0

is also called a sigmoid function. The graph of u(t). is obtained from the graph of
S(t). by a sequence of three linear transformations: (i) shift to the left by α0 ., (ii) scale
by a diagonal matrix diag{1/α1 , β1 }., and (iii) shift up by β0 .. In this book, we shall
scale S(t). by a diagonal matrix diag{1/2, 2}. and then shift it down by 1; namely,
we choose α0 = 0., α1 = β1 = 2., and β0 = −1.. The resulting function is the same
as the hyperbolic tangent function:

2 1 − e−2t et − e−t
u(t) = 2S(2t) − 1 =
. − 1 = = = tanh t. (4.5)
1 + e−2t 1 + e−2t et + e−t

It is easily seen that u(t) = tanh t . is an odd and increasing function bounded by two
horizontal asymptotes ± 1.. We also note that

1
u (t) =
. = 1 − tanh2 t = 1 − [u(t)]2 . (4.6)
cosh2 t

The sigmoid function S(t). and its linear transformations are illustrated in Fig. 4.1.

2
S(t)
2S(2t)
1.5 tanh(t)

0.5

-0.5

-1
-5 0 5

Fig. 4.1 The sigmoid function S(t) = 1/(1 + e−t ). and its linear transformations
48 4 Nonlinear Regression

4.3 Optimization Problem

In the nonlinear regression model, the prediction vector is a nonlinear function

Y = tanh(xp),
. (4.7)

where, as usual, x ∈ Rn×(d+1) . is the input matrix, and p ∈ Rd+1 . is the model
parameter to be estimated. The above formula is interpreted as Yi = tanh((xp)i ). for
i = 1, · · · , n.. The output vector is binary: y ∈ {−1, 1}n ., which is different from
that in the linear regression model. Another difference between nonlinear and linear
regression model lies in the choice of objective function. In the linear regression
model, the Euclidean distance is used to measure the error between the prediction
vector and the observation vector. In the nonlinear regression model, we increase the
penalty when the prediction is closer to the opposite side of the observation. Given
a binary observation v ∈ {−1, 1}. and a prediction u ∈ (−1, 1)., we define the cost
function
2 ln 2 − (1 + v) ln(1 + u) − (1 − v) ln(1 − u)
C(u, v) :=
. . (4.8)
2
It is easily seen that

∂ 1+v 1−v u−v

. C(u, v) = − + = . (4.9)
∂u 2(1 + u) 2(1 − u) 1 − u2

Note that the value of v is binary: either v = 1. or v = −1.. An equivalent definition

of the cost function is

2 ln[2/(1 + u)], v = 1,
C(u, v) = ln
. = (4.10)
1 + uv ln[2/(1 − u)], v = −1.

We consider the following optimization problem.

Optimization Problem
Let x ∈ Rn×(d+1) . and y ∈ {−1, 1}n . be given. Find p ∈ Rd+1 . such that the
objective function
n
J (p) :=
. C(Yi , yi ) (4.11)
i=1

is minimized. Here, Yi = tanh((xp)i ). for i = 1, · · · , n., and C is the cost

function defined in (4.8).
4.4 Gradient Descent Method 49

4.4 Gradient Descent Method

Let J (p). be the objective function given in (4.11). We will use the gradient descent
method to find a numerical approximation of the parameter vector p that minimizes
J (p).. For any i = 1, · · · , n., we obtain from
⎛ ⎞
d
Yi = tanh((xp)i ) = tanh ⎝
. xij pj ⎠
j =0

that
∂Yi
. = (1 − Yi2 )xij , j = 0, · · · , d.
∂pj

It then follows from (4.9) and (4.11) that

n
∂
. J (p) = (Yi − yi )xij ,
∂pj
i=1

for j = 0, · · · , d .. The above equation can be written in matrix form as

J (p) = x T (Y − y).
. (4.12)

The following lemma motivates the gradient descent method.

Lemma 4.1 Given p ∈ Rd+1 . and J ∈ C 2 (U, R)., where U ⊂ Rd+1 . is an open
neighborhood of p. Assume J (p) 0.. Then there exists α0 > 0. (which generally
depends on p) such that J (p − αJ (p)) < J (p). for any α ∈ (0, α0 )..
Proof Since J ∈ C 2 (U, R). and J (p) 0., we can find a small α0 > 0. such that
p − αJ (p) ∈ U . and

α|[J (p)]T J (p − αJ (p))J (p)| < J (p)

.
2
2

for all α ∈ (0, α0 ).. It then follows from Taylor expansion that

α2
J (p − αJ (p)) − J (p) = −α J (p)
.
2
2 + [J (p)]T J (p − ξ J (p))J (p),
2

for some ξ ∈ (0, α).. Coupling the above two inequalities gives J (p − αJ (p)) <
J (p).. This completes the proof.
The gradient descent method is an iterative method to find the minimizer along the
direction of descent gradient; namely,
50 4 Nonlinear Regression

p(k+1) = p(k) − αJ (p(k) ) = p (k) − αx T [tanh(xp (k) ) − y],

. (4.13)

where α > 0. is called the learning rate, which is a hyperparameter to be tuned in

practical applications.

4.5 Python Code

We use the nonlinear regression model to fit the AND gate. The gradient descent
method is implemented in the following Python code:

Nonlinear Regression Code

import numpy as np

def nonlinear_regression(x,y,al,iter_max):
d=[Link][1]-1
p=[Link](0,1,(3,1))
for iter in range(iter_max):
p=p-al*[Link](x.T,[Link]([Link](x,p))-y)
return p

x=[Link]([[1,1,1],[1,1,-1],[1,-1,1],[1,-1,-1]])
y=[Link]([[1],[-1],[-1],[-1]])

p=nonlinear_regression(x,y,al=0.1,iter_max=1000)
print(p)
Y=[Link]([Link](x,p))
print([Link]((Y,y),axis=1))

The result of numerical iteration is given as below:

Iteration Result
[[-2.99576037]
[ 2.99576042]
[ 2.99576041]]
[[ 0.99501275 1. ]
[-0.99501275 -1. ]

(continued)
4.6 Discussions 51

[-0.99501275 -1. ]
[-0.99999997 -1. ]]

It is noted that the nonlinear regression model fits the AND gate with a high
accuracy.

4.6 Discussions

We first discuss how the nonlinear regression model works for the AND gate. Note
from the result of numerical iteration via the gradient descent method in Sect. 4.5
that p = (−M, M, M)T ∈ R3 ., where M ≈ 3.. A simple calculation gives
⎛ ⎞ ⎛ ⎞
1 1 1 ⎛ ⎞ M
⎜1 1 −1⎟ −M ⎜ ⎟
.xp = ⎜ ⎟ ⎝ M ⎠ = ⎜ −M ⎟ .
⎝1 −1 1 ⎠ ⎝ −M ⎠
M
1 −1 −1 −3M

Hence, the prediction vector is

⎛ ⎞
tanh(M)
⎜ − tanh(M) ⎟
.Y = tanh(xp) = ⎜ ⎟
⎝ − tanh(M) ⎠ .
− tanh(3M)

Note from Fig. 4.1 that the function tanh(t). has two horizontal asymptotes ± 1.;
namely, tanh(t) → ±1. as t → ±∞.. This implies that tanh(M) ≈ 1. for large
M > 0.. Actually, when M = 3., we have tanh(M) ≈ 0.995. and tanh(3M) ≈
0.99999997.. Hence, the prediction vector Y calculated in the above formula can be
regarded as a good approximation of the observation vector y = (1, −1, −1, −1)T ..
Finding the global minimum of a general objective function, with or without
constraints, remains an open problem in optimization. The effectiveness of machine
learning lies in the fact that an accurate approximation of the parameter vector need
not correspond to a global minimum, which in some cases may not even exist.
Instead, there often exist infinitely many parameter vectors that yield sufficiently
good approximations, and machine learning algorithms are able to identify such
solutions efficiently.
For example, in the nonlinear regression model of the AND gate, the parameter
vector p = (−M, M, M)T . is a good solution for any sufficiently large M > 0.. This
illustrates that, in practice, it is neither necessary nor typically possible to determine
the global minimum in a machine learning model. Rather, it suffices to perform
52 4 Nonlinear Regression

a finite number of numerical iterations, such as those generated by the gradient

descent method, to obtain predictions that closely approximate the observed data.
In the case of the nonlinear regression model for the AND gate, we show in the
following proposition that the global minimum of the objective function J (p). in
(4.11) cannot be attained at any finite parameter vector.
Proposition 4.1 Consider the AND gate with input matrix x and output vector y
given as in (4.1) and (4.2), respectively. Let J (p). be the objective function defined
in (4.11). We have J (p) > 0. for all p ∈ R3 .. Moreover, we have J (q(M)) → 0. as
M → ∞., where q(M) := (−M, M, M)T ∈ R3 ..
Proof Note from (4.5) that C(u, v) > 0. for any u ∈ (−1, 1). and v ∈ {−1, 1}..
Hence, we obtain from (4.11) that J (p) > 0. for all p ∈ R3 .. Moreover, a simple
calculation gives
⎛ ⎞ ⎛ ⎞
tanh(M) 1
⎜ − tanh(M) ⎟ ⎜−1⎟
. tanh(xq(M)) = ⎜ ⎟ ⎜ ⎟
⎝ − tanh(M) ⎠ → ⎝−1⎠ = y
− tanh(3M) −1

as M → ∞.. Thus,

4
. lim J (q(M)) = C(yi , yi ) = 0.
M→∞
i=1

This completes the proof.

The above proposition implies that the optimization problem for the nonlinear
regression model of the AND gate is not well-posed because the global minimum
does not exist. The amazing thing in machine learning is that the gradient descent
method can still find a good approximation efficiently even though the global
minimum may not exist. The open problem of global optimization is not solved
but effectively avoided in machine learning.
Nonlinear regression and linear regression differ in three aspects: the observation
vector, the model, and the cost function. However, the derivatives of the objective
function and the iterative functions of the gradient descent method for both
regression models have the similar formulas. A comparison of these two models
is listed in Table 4.2.

Problems
4.1 Use a linear regression model to fit the OR gate. Find the parameter vector and
the prediction vector.
4.2 Verify that the function u(t) = (1 + e−2t )−1 . satisfies the differential equation
u (t) = u(t)[1 − u(t)]. with initial condition u(0) = 1/2..
4.6 Discussions 53

Table 4.2 The comparison between nonlinear regression and linear regression
Linear regression Nonlinear regression
Input matrix x ∈ Rn×(d+1) . x ∈ Rn×(d+1) .
Observation vector y ∈ Rn . y ∈ {−1, 1}n .
Model Y = xp ∈ Rn . Y = tanh(xp) ∈ (−1, 1)n .
Parameter vector p ∈ Rd+1 . p ∈ Rd+1 .
Cost function C(u, v) = (u − v)2 /2. C(u, v) = ln[2/(1 + uv)].
Objective function J (p) = ni=1 C(Yi , yi ). J (p) = ni=1 C(Yi , yi ).
Gradient J (p) = x T (Y − y). J (p) = x T (Y − y).
Gradient descent method p (k+1) = p(k) − αJ (p (k) ). p (k+1) = p(k) − αJ (p (k) ).

4.3 Verify that the function u(t) = (1 + e−2t )−1 . is the cumulative distribution
function of the logistic distribution with location parameter 0 and scale parameter
1/2.; namely,
t
u(t) =
. v(s)ds,
−∞

where v(s) = 2e−2s (1 + e−2s )−2 . is the probability density function.

4.4 List the sequence of linear transformations from u(t) = tanh t . to S(t) = (1 +
e−t )−1 ..
4.5 Let v ∈ {−1, 1}., u ∈ (−1, 1)., and C(u, v). be defined in (4.8). Prove (4.10).
4.6 Let J (p). be the objective function given in (4.11). Verify that

J (p) = x T diag(1 − Y 2 )x + λI,

.

where I is the identity matrix in R(d+1)×(d+1) . and diag(1 − Y 2 ) ∈ Rn×n . is a

diagonal matrix with diagonal terms being 1 − Yi2 . for i = 1, · · · , n.. Moreover,
prove that J (p). is positive definite.
4.7 Given p ∈ Rd+1 . and J ∈ C 2 (U, R)., where U ⊂ Rd+1 . is an open
neighborhood of p, assume J (p) 0., and choose a small α0 > 0. such that U
contains the open ball Bα0 J (p) 2 (p). centered at p with radius α0 J (p) 2 .. For any
q = p − αJ (p). with α ∈ (0, α0 )., prove that

α2
J (q) − J (p) = −α J (p)
.
2
2 + [J (p)]T J (p − ξ J (p))J (p),
2
for some ξ ∈ (0, α)..
4.8 Use the nonlinear regression model to fit the OR gate; namely, find the
numerical approximation of p that minimizes the objective function J (p). in (4.11),
54 4 Nonlinear Regression

where the input matrix x is given as in (4.1), and the output vector is y =
(1, 1, 1, −1)T ∈ R4 ..
4.9 Consider the model Y = tanh(xp)., where x is given as in (4.1). Let M > 0. be a
sufficiently large constant. Determine which logic gate does the model approximate
when p is chosen as each of the following seven vectors:
⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞
M M M M −M −M −M
. ⎝M ⎠ , ⎝ M ⎠ , ⎝−M ⎠ , ⎝−M ⎠ , ⎝ M ⎠ , ⎝−M ⎠ , ⎝−M ⎠ .

M −M M −M −M M −M

4.10 Given J ∈ C 2 (Rd+1 , R). such that J (p). is positive definite for all p ∈ Rd+1 .,
prove that p ∈ Rd+1 . is a global minimum of J if and only if it is a critical point of
J ; namely, J (p) < J (q). for all q ∈ Rd+1 . if and only if J (p) = 0..
Chapter 5
Shallow Neural Network

Abstract The nonlinear model introduced in the previous chapter consists of

applying a sigmoid function to a linear combination of the input data. This model
can successfully approximate several logic gates, including the AND, OR, NAND,
and NOR gates. However, it fails to represent more complex gates such as the XOR
gate (as well as the XAND gate). In this chapter, we introduce intermediate variables
and employ composite functions to approximate the XOR gate. These composite
functions can be interpreted as a network of interconnected neural nodes, and the
intermediate variables correspond to the hidden layers of the network.

5.1 Notations

For convenience, we introduce some notations to denote the operations of logic

gates. As in the previous chapter, we use 1 and − 1. to denote True and False,
respectively. We are considering operations on a binary space B = {1, −1}.. For any
x ∈ B ., the negative (opposite) of x is denoted by − x .; namely, False is the negative
of True and vice versa.
Let x1 ∈ B . and x2 ∈ B . be two input variables of a logic gate. The AND gate is
the same as the minimum of x1 . and x2 ., which is denoted as

1, x1 = x2 = 1,
x1 ∧ x2 = min{x1 , x2 } =
. (5.1)
−1, (x1 + 1)(x2 + 1) = 0.

The OR gate is the same as the maximum of x1 . and x2 ., which is denoted as

1, (x1 − 1)(x2 − 1) = 0,
.x1 ∨ x2 = max{x1 , x2 } = (5.2)
−1, x1 = x2 = −1.

The XAND gate is the same as the product of x1 . and x2 ., which is denoted as

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 55

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
56 5 Shallow Neural Network

1, x1 = x2 ,
x1 · x2 =
. (5.3)
−1, x1 x2 .

The XOR gate is the negative of the XAND gate: − (x1 · x2 ).. One can easily verify
the following identities:

. − (x1 ∧ x2 ) =(−x1 ) ∨ (−x2 ), . (5.4)

− (x1 ∨ x2 ) =(−x1 ) ∧ (−x2 ), . (5.5)
− (x1 · x2 ) =(−x1 ) · x2 = x1 · (−x2 ). (5.6)

5.2 The XOR Gate

A logic gate has two binary inputs x1 . and x2 . and one binary output y. The
corresponding data has a dimension of d = 2. and a size of n = 4.. Table 5.1
illustrates the input and output data of the AND gate, OR gate, NAND gate, NOR
gate, and XOR gate, respectively.
The input matrix and the output vector of the XOR gate can be written as
⎛ ⎞
1 1 1
⎜1 1 −1⎟
.x = ⎜ ⎟
⎝1 −1 1 ⎠ ∈ R ,
4×3
(5.7)
1 −1 −1

and
⎛ ⎞
−1
⎜1⎟
.y = ⎜ ⎟
⎝ 1 ⎠∈R ,
4
(5.8)
−1

respectively. If we run the nonlinear regression code in Sect. 4.5 with the above x
and y, the simulation result is given below:

Table 5.1 The logic gates x1 . x2 . AND OR NAND NOR XOR

with 1 and − 1. denoting True
and False, respectively 1 1 1 1 −1
. −1
. −1
.

1 −1. −1
. 1 1 −1
. 1
−1. 1 −1
. 1 1 −1
. 1
−1. −1. −1
. −1
. 1 1 −1
.
5.3 Composite Function 57

Simulation Result
[[ 1.08017916e-17]
[ 7.93636823e-18]
[-2.85865445e-17]]
[[-9.84838471e-18 -1.00000000e+00]
[ 4.73247043e-17 1.00000000e+00]
[-2.57211212e-17 1.00000000e+00]
[ 3.14519678e-17 -1.00000000e+00]]

The parameter vector is close to (0, 0, 0)T ., and the prediction vector is also close
to the zero vector (0, 0, 0, 0)T ., which fails to approximate the observation vector of
XOR gate given in (5.8).

5.3 Composite Function

A careful examination of Table 5.1 shows that the XOR gate can be regarded as the
AND gate whose inputs are the outputs of the NAND gate and the OR gate; namely,
for any x1 ∈ B . and x2 ∈ B ., we have

. − x1 · x2 = [−(x1 ∧ x2 )] ∧ (x1 ∨ x2 ). (5.9)

Since each of the AND, NAND, and OR gates can be approximated by a nonlinear
sigmoid function with a linear combination, we approximate y as a composite
function of x1 . and x2 . in the sense that y is a function of two intermediate variables
z1 . and z2 . that approximates the AND gate, where z1 . and z2 . are functions of x1 . and
x2 . that approximate the NAND gate and the OR gate, respectively. This suggests us
to propose the following model. For each i = 1, 2, 3, 4., we define

Yi = tanh(zi0 p0 + zi1 p1 + zi2 p2 ),

. (5.10)

where zi0 = 1. and

.zi1 = tanh(xi0 q0 + xi1 q1 + xi2 q2 ), . (5.11)

zi2 = tanh(xi0 r0 + xi1 r1 + xi2 r2 ). (5.12)

Note that the composite function has a total of nine parameters: pj ., qj ., and rj . with
j = 0, 1, 2.. The objective function is given by
58 5 Shallow Neural Network

n
J (p, q, r) :=
. C(Yi , yi ), (5.13)
i=1

where C is the cost function defined in (4.8). We also recall from (4.9) that

∂C Yi − yi
. = .
∂Yi 1 − Yi2

We need to calculate the partial derivatives Jp , Jq , Jr ∈ R3 . so as to implement the

gradient descent method. By the chain rule, we obtain the j -th component of Jp . as

n n
∂J ∂C ∂Yi
. = = (Yi − yi )zij , (5.14)
∂pj ∂Yi ∂pj
i=1 i=1

for each j = 0, 1, 2.. Similarly, we have the j -th components of Jq . and Jr .,

respectively,
n n
∂J ∂C ∂Yi ∂zi1
. = = (Yi − yi )p1 (1 − zi1
2
)xij , . (5.15)
∂qj ∂Yi ∂zi1 ∂qj
i=1 i=1
n n
∂J ∂C ∂Yi ∂zi2
= = (Yi − yi )p2 (1 − zi2
2
)xij , (5.16)
∂rj ∂Yi ∂zi2 ∂rj
i=1 i=1

for j = 0, 1, 2.. In matrix form, we denote the gradient vectors as

∂J ∂J ∂J T ∂J ∂J ∂J T ∂J ∂J ∂J T
Jp = (
. , , ) , Jq = ( , , ) , Jr = ( , , )
∂p0 ∂p1 ∂p2 ∂q0 ∂q1 ∂q2 ∂r0 ∂r1 ∂r2

and obtain

Jp =zT (Y − y), .
. (5.17)
Jq =x T [(Y − y) ∗ (1 − z1 ∗ z1 )]p1 , . (5.18)
Jr =x T [(Y − y) ∗ (1 − z2 ∗ z2 )]p2 , (5.19)

where z1 = (z11 , · · · , z41 )T ., z2 = (z12 , · · · , z42 )T ., and ∗. denotes the component-

wise multiplication; namely, if A ∈ Rn×m . and B ∈ R n×m ., then A ∗ B ∈ Rn×m . such
that (A ∗ B)ij = Aij Bij . for each i = 1, · · · , n. and j = 1, · · · , m.. We then write
the gradient descent method as

p(k+1) = p(k) − αJp(k) , q (k+1) = q (k) − αJq(k) , r (k+1) = r (k) − αJr(k) ,

.

(5.20)

where α > 0. is the learning rate.

5.5 Discussions 59

5.4 Python Code

The Python code for the gradient descent method is given below:

XOR Code
import numpy as np

x=[Link]([[1,1,1],[1,1,-1],[1,-1,1],[1,-1,-1]])
y=[Link]([[-1],[1],[1],[-1]])

al=0.1
iter_max=1000
p=[Link](0,1,(3,1))
q=[Link](0,1,(3,1))
r=[Link](0,1,(3,1))

for iter in range(iter_max):

z1=[Link]([Link](x,q))
z2=[Link]([Link](x,r))
z=[Link](([Link]((4,1)),z1,z2),axis=1)
Y=[Link]([Link](z,p))
Jp=[Link](z.T,Y-y)
Jq=[Link](x.T,(Y-y)*(1-z1*z1))*p[1,0]
Jr=[Link](x.T,(Y-y)*(1-z2*z2))*p[2,0]
p=p-al*Jp
q=q-al*Jq
r=r-al*Jr

print([Link]((p,q,r),axis=1))
print([Link]((Y,y),axis=1))

5.5 Discussions

A sample of simulation result is given as below:

60 5 Shallow Neural Network

Simulation Result
[[-3.08739711 1.80172906 1.5999076 ]
[ 3.4102806 -1.95476405 1.77937427]
[ 3.42289907 -1.94572899 1.77547159]]
[[-0.99478607 -1. ]
[ 0.99725619 1. ]
[ 0.99726898 1. ]
[-0.99469775 -1. ]]

From the simulation result, we have the following observations:

• p ≈ (−3, 3, 3)T ., which is an approximation of the AND gate: Y ≈ z1 ∧ z2 ..
• q ≈ (2, −2, −2)T ., which is an approximation of the NAND gate: z1 ≈ −(x1 ∧
x2 )..
• r ≈ (2, 2, 2)T ., which is an approximation of the OR gate: z2 ≈ x1 ∨ x2 ..
This simulation result agrees with the identity (5.9). However, since the decompo-
sition of the XOR gate is not unique, the parameter values obtained by gradient
descent method with a random initial guess may differ. For instance, another
simulation gives p ≈ (3, 3, 3)T ., q ≈ (−2, 2, −2)T ., and r ≈ (−2, −2, 2)T ., which
approximates the following decomposition:

. − x1 · x2 = [x1 ∧ (−x2 )] ∨ [(−x1 ) ∧ x2 ]. (5.21)

This example demonstrates that the optimization problem may lack a finite exact
solution and can have multiple approximate solutions. Additionally, the objective
function is non-convex. It is important to note that, depending on the initial guess,
the gradient descent method may sometimes converge to a critical point that does
not minimize the objective function.
Our model (5.10)–(5.12) can be interpreted as a shallow neural network x →
z → Y ., where x ∈ R4×3 . is the input layer, Y ∈ R4 . is the output layer, and z ∈ R4×3 .
is the hidden layer. The map from x to z is characterized by an activation (hyperbolic
tangent) function together with linear combinations of x determined by a parameter
matrix (q, r) ∈ R3×2 ., while the map from z to Y is characterized by an activation
(hyperbolic tangent) function together with a linear combination of z determined by
a parameter vector p ∈ R3 ..

Problems
5.1 Given x1 , x2 ∈ B = {1, −1}., prove that x1 ∨x2 = x2 ∨x1 . and x1 ∧x2 = x2 ∧x1 ..
5.2 Given x1 , x2 , x3 ∈ B = {1, −1}., prove that (x1 ∨ x2 ) ∨ x3 = x1 ∨ (x2 ∨ x3 ).
and (x1 ∧ x2 ) ∧ x3 = x1 ∧ (x2 ∧ x3 )..
5.3 Given x1 , x2 , x3 ∈ B = {1, −1}., prove that
5.5 Discussions 61

(x1 ∨ x2 ) ∧ x3 = (x1 ∧ x3 ) ∨ (x2 ∧ x3 ),

.

and

. (x1 ∧ x2 ) ∨ x3 = (x1 ∨ x3 ) ∧ (x2 ∨ x3 ).

5.4 Verify (5.9) and (5.21).

5.5 Prove that the zero vectors p = (0, 0, 0)T ., q = (0, 0, 0)T ., and r = (0, 0, 0)T .
remain invariant under the gradient descent iteration (5.20).
5.6 Given

12
A=
. ,
34

compute A2 . and A ∗ A..

Chapter 6
Deep Neural Network

Abstract We generalize the shallow neural network introduced in the previous

chapter to a deep neural network. In mathematical terms, a deep neural network can
be viewed as a nested composite function comprising multiple layers of intermediate
(matrix-valued) variables and several parameter matrices that define the functions
connecting consecutive layers. The gradient of the objective function with respect to
each parameter matrix is obtained through a straightforward application of the chain
rule in multivariate calculus. This process, which proceeds in the reverse order of the
forward computation of the objective function, is known as backward propagation.

6.1 The Model

Recall that the model of shallow neural network (5.10)–(5.12) is rewritten as a chain
of maps x → z → Y ., where x ∈ R4×3 . is the input layer, Y ∈ R4 . is the output layer,
and z ∈ R4×3 . is the hidden layer. Now we use a unified notation to denote a deep
neural network

x (0) → x (1) → · · · → x (L−1) → x (L) ,

. (6.1)

where x (0) = x ∈ Rn×(d+1) . is the input layer and x (L) = (1n , Y ) ∈ Rn×(1+1) . with
1n := (1, · · · , 1)T ∈ Rn . is the output layer, while for each l = 1, · · · , L − 1.,
x (l) ∈ Rn×(dl +1) . with x10 = · · · = xn0 = 1. is the hidden layer. The dimension
(l) (l)

of the l-th layer (with l = 0, 1, · · · , L.) is denoted by dl .; in particular, d0 = d .

is the number of input features, and dL = 1. is the dimension of the output layer.
The depth of the neural network is denoted by L. The map x (l) → x (l+1) . (with
l = 0, 1, · · · , L − 1.) is defined as

x (l+1) = (1n , tanh(x (l) p(l) )),

. (6.2)

where p(l) ∈ R(dl +1)×dl+1 . is the parameter matrix connecting the l -th layer with the
(l + 1).-th layer. In other words, the above map can be written as

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 63

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
64 6 Deep Neural Network

(l+1) 1, jl+1 = 0,
xijl+1 =
.
dl (l) (l) (6.3)
tanh( jl =0 xijl pjl jl+1 ), jl+1 = 1, · · · , dl+1 ,

for each i = 1, · · · , n., l = 0, · · · , L − 1., and jl+1 = 0, 1, · · · , dl+1 .. The objective

function is given by
n
J (p(0) , · · · , p (L−1) ) :=
. C(Yi , yi ), (6.4)
i=1

where C is the cost function defined in (4.8).

6.2 Backward Propagation

Similar to the case of shallow neural network, we shall calculate the gradients
∂J /∂p(l) . for each l = 0, · · · , L − 1.. For convenience, we introduce the following
notations:

u(l) :=x (l) p(l) ∈ Rn×dl+1 , .

. (6.5)
∂J
v (l) := ∈ Rn×dl+1 , . (6.6)
∂u(l)
∂J
q (l) := (l) ∈ R(dl +1)×dl+1 , (6.7)
∂p

for each l = 0, · · · , L − 1.. The following formulas are elementary but crucial:

(l) (l) (l)

dxijl (l)
∂uijl+1 (l)
∂uijl+1 (l)
.
(l−1)
= 1 − [xijl ]2 , (l)
= pjl jl+1 , (l)
= xijl , (6.8)
duijl ∂xijl ∂pjl jl+1

for each i = 1, · · · , n., l = 0, · · · , L − 1., jl = 0, · · · , dl ., and jl+1 = 1, · · · , dl+1 ..

Now, we are ready to state the backward propagation:

(L) (L)
(L−1) ∂J ∂xi1 xi1 − yi (L)
vi1
. = (L) (L−1)
= (L)
{1 − [xi1 ]2 } = Yi − yi , . (6.9)
∂xi1 ∂ui1 1 − [xi1 ]2
(l)
(l−1)
dl+1
∂J ∂uijl+1 dxij(l)l dl+1
(l ) (l) (l)
vijl = (l) (l) (l−1)
= vijl+1 pjl jl+1 {1 − [xijl ]2 },
jl+1 =1 ∂uijl+1 ∂xijl duijl jl+1 =1
(6.10)
6.3 Python Code 65

Table 6.1 Neural network notations

Notation Dimension Remark Relation
x (l) . n × (dl + 1). Layer x (l+1) = (1n , tanh(u(l) )).
p (l) . (dl + 1) × dl+1 . Parameter matrix u(l) = x (l) p (l) .
u(l) . n × dl+1 . Linear combination u(l) = x (l) p (l) .
v (l) . n × dl+1 . Gradient v (l) = ∂J /∂u(l) .
q (l) . (dl + 1) × dl+1 . Gradient q (l) = ∂J /∂p(l) = [x (l) ]T v (l) .

for each i = 1, · · · , n., l = L − 1, · · · , 1., and jl = 1, · · · , dl .. Note that the last

formula provides a backward recurrence relation for v (l) .. Finally, we have

n (l) n
(l) ∂J ∂uijl+1 (l) (l)
.q
jl jl+1 = (l) (l)
= vijl+1 xijl , (6.11)
i=1 uijl+1 ∂pjl jl+1 i=1

for each l = 0, · · · , L − 1., jl = 0, · · · , dl ., and jl+1 = 1, · · · , dl+1 .. To write the

above formulas in matrix form, we shall introduce several notations. Denote by ∗.
the componentwise matrix multiplication. Let 1n×(dl +1) ∈ Rn×(dl +1) . be the matrix
with ones on all entries. For any A ∈ Rn×(dl +1) ., we denote by Ā(1:n)×(1:dl ) ∈ Rn×dl .
the submatrix of A obtained by removing the column jl = 0.. We can rewrite the
above formulas as

v (L−1) =Y − y, .
. (6.12)
(1:n)×(1:dl )
v (l−1) ={v (l) [p (l) ]T } ∗ [1n×(dl +1) − x (l) ∗ x (l) ] , l = L − 1, · · · , 1, .
(6.13)
q (l) =[x (l) ]T v (l) , l = 0, · · · , L − 1. (6.14)

The gradient descent method is thus to replace p (l) . with p(l) − αq (l) . during each
iteration, where α > 0., as usual, denotes the learning rate. The neural network
notations are illustrated in Table 6.1.

6.3 Python Code

The Python code for the (deep) neural network is given below:
66 6 Deep Neural Network

Neural Network Code

import numpy as np

def neural_network(x,y,hidden_layers,al,iter_max):
n,m=[Link]
d=[m-1]+hidden_layers+[1]
L=len(d)-1
x=[x]+[None]*L
v=[None]*L
p=[None]*L
for l in range(L):
p[l]=[Link](0,1,(d[l]+1,d[l+1]))
for iter in range(iter_max):
for l in range(L):
u=[Link]([Link](x[l],p[l]))
x[l+1]=[Link](([Link]((n,1)),u),axis=1)
Y=x[L][:,1].reshape(n,1)
v[L-1]=Y-y
for l in range(L-1,0,-1):
v[l-1]=([Link](v[l],p[l].T)*(1-x[l]*x[l]))[:,1:]
for l in range(L):
p[l]=p[l]-al*[Link](x[l].T,v[l])
return Y,p

x=[Link]([[1,1,1],[1,1,-1],[1,-1,1],[1,-1,-1]])
y=[Link]([[-1],[1],[1],[-1]])
Y,p=neural_network(x,y,hidden_layers=[2],al=0.1,
iter_max=1000)
print(p)
print([Link]((Y,y),axis=1))

6.4 Discussions

In the Python code, we use the list data type to address the layers.
1. d is a list of L + 1. numbers, where d[l]. with l = 0, · · · , L. corresponds to the
dimension in l-th layer. In particular, d[0]. is the dimension of the input layer, and
d[L] = 1. is the dimension of the output layer.
6.4 Discussions 67

2. x is a list of L + 1. matrices, where x[l] ∈ Rn×(dl +1) . with l = 0, · · · , L.

corresponds to the data of l-th layer. In particular, x[0]. is the input data, and
y[0] = (1n , Y ). is the prediction on the output layer.
3. p is a list of L matrices, where p[l] ∈ R(dl +1)×dl+1 . with l = 0, · · · , L − 1.
corresponds to the parameter matrix connecting the l-th layer to the l + 1.-th
layer.
4. v is a list of L matrices, where v[l] ∈ Rn×dl+1 . with l = 0, · · · , L−1. corresponds
to the gradient of J with respect to u.
The four major steps in the main function are:
1. Initial random guess of parameter matrices:

p[l] = [Link](0, 1, (d[l] + 1, d[l + 1]))

.

2. Forward propagation:

x[l + 1] = [Link](([Link]((n, 1)), u), axis = 1)

.

3. Backward propagation:

v[l − 1] = ([Link] (v[l], p[l].T ) ∗ (1 − x[l] ∗ x[l]))[:, 1 :]

.

4. Gradient descent iteration:

.p[l] = p[l] − al ∗ [Link] (x[l].T , v[l])

A sample of simulation result is given as below:

Simulation Result
[array([[-1.5815611 , 1.82134624],
[ 1.76355366, 1.97592353],
[-1.75918357, -1.96452545]]), array([[ 3.04703861],
[ 3.38373003],
[-3.36820299]])]
[[-0.997053 -1. ]
[ 0.99426201 1. ]
[ 0.99437418 1. ]
[-0.9970351 -1. ]]

With only one hidden layer of dimension 2, the simulation results align with
those of the shallow neural network discussed in the previous chapter. To adapt this
68 6 Deep Neural Network

code for a deeper neural network, we simply need to modify the list hidden_layer
by increasing the number of hidden layers and configuring them with varying
dimensions.

Problems
6.1 Let C = AB ., where A ∈ Rn×m ., B ∈ Rm×r ., and C ∈ Rn×r .. For each i =
1, · · · , n., j = 1, · · · , r ., k = 1, · · · , n., l = 1, · · · , m., s = 1, · · · , r ., prove that

∂Cij ∂Cij
. = δik Blj , = δj s Ail ,
∂Akl ∂Bls

where δik . (resp. δj s .) is the Kronecker delta symbol which equals 1 when i = k .
(respectively, j = s .) and 0 when i k . (respectively, j s .).
6.2 Verify (6.8).
6.3 Verify (6.9) and (6.10).
6.4 Verify (6.13).
6.5 Change the hidden layers in the Python code, and observe the simulation
results.
Chapter 7
Batch Normalization

Abstract This chapter introduces batch normalization in deep neural networks, a

technique designed to stabilize and accelerate training by normalizing the interme-
diate variables at each layer. We extend the deep neural network model from the
previous chapter, where each layer consists of a linear combination of the previous
layer’s outputs followed by an activation function, by adding a normalization step
prior to the linear combination. Specifically, for each hidden layer, we standardize
the activations using the layer-wise mean and standard deviation. The chapter
also presents the backward propagation algorithm adapted to batch normalization,
detailing how gradients propagate through the normalization step. Key formulas
for derivatives with respect to parameters, activations, and normalized variables are
derived, highlighting the dependence of each output on all samples in the batch.

7.1 The Model

As in the previous chapter, we consider a neural network denoted by

x (0) → x (1) → · · · → x (L−1) → x (L) ,

. (7.1)

where x (0) = x ∈ Rn×(d0 +1) . is the input layer and x (L) = (1n , Y ) ∈ Rn×(dL +1) .
(with dL = 1. and 1n = (1, · · · , 1)T ∈ Rn .) is the output layer, while for each
l = 1, · · · , L − 1., x (l) ∈ Rn×(dl +1) . is the hidden layer. The dimension of the l-th
layer (with l = 0, 1, · · · , L.) is denoted by dl .; in particular, d0 . is the number of
input features, and dL = 1. is the dimension of the output layer. The depth of the
neural network is denoted by L. The map x (l) → x (l+1) . (with l = 0, 1, · · · , L − 1.)
was defined in the previous chapter as a combination of an activation (generalized
sigmoid/hyperbolic tangent) function and a linear combination of x (l) . involving a
parameter matrix p(l) ∈ Rdl +1×dl+1 .. The key idea of batch normalization is to
introduce an additional step of normalization before taking the linear combination.
In other words, for each l = 0, 1, · · · , L − 1., i = 1, · · · , n., and jl = 0, · · · , dl ., we
define

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 69

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
70 7 Batch Normalization

(l) 1, jl = 0,
zijl =
.
(l) (l) (l) (7.2)
[xijl − μjl ]/σjl , jl = 1, · · · , dl ,

where
n
(l) 1 (l)
μjl =
. xijl (7.3)
n
i=1

is the average and

n
(l) 1 (l) (l)
σj l =
. [xijl − μjl ]2 (7.4)
n
i=1

is the standard deviation. Next, we set

u(l) =z(l) p(l) ∈ Rn×dl+1 , .

. (7.5)
x (l+1) =(1n , tanh(u(l) )), (7.6)

where p(l) ∈ R(dl +1)×dl+1 . is the parameter matrix connecting the l -th layer with the
(l + 1).-th layer, and 1n ∈ Rn . is the vector of 1s. The objective function is still given
by
n
J (p(0) , · · · , p (L−1) ) :=
. C(Yi , yi ), (7.7)
i=1

where C is the cost function defined in (4.8).

7.2 Backward Propagation

For convenience, we introduce the following notations:

∂J
v (l) :=
. ∈ Rn×dl+1 , . (7.8)
∂u(l)
∂J
w (l) := (l) ∈ Rn×(dl +1) , . (7.9)
∂z
∂J
q (l) := (l) ∈ R(dl +1)×dl+1 , (7.10)
∂p
7.2 Backward Propagation 71

for each l = 0, · · · , L − 1.. The major difference introduced by normalization is that

u(l)
ijl+1 . not only depends on uijl
(l−1)
. but also depends on u
(l−1)
i jl . for all i = 1, · · · , n..
For each l = 0, · · · , L − 1., i = 1, · · · , n., i = 1, · · · , n., and jl = 1, · · · , dl ., we
obtain from (7.2), (7.3), and (7.4) that

∂μ(l)
jl
(l)
1 ∂σjl xij(l)l − μ(l)
jl
(l)
zij
.
(l)
= , = = l
, (7.11)
∂xijl n ∂x (l) nσ
(l) n
ijl jl

and
(l) (l) (l) (l) (l) (l)
∂zi jl δii − 1/n xi jl − μjl zijl nδii − 1 − zijl zi jl
.
(l)
= (l)
− (l)
· = (l)
, (7.12)
∂xijl σj l [σjl ]2 n nσjl

where δii . is the Kronecker delta symbol which equals 1 when i = i . and 0 when
i i ..
Similar to (6.8), we also need the following formulas:

(l) (l)
∂uijl+1 (l)
∂uijl+1 (l)
.
(l)
= pjl jl+1 , (l)
= zijl , (7.13)
∂zijl ∂pjl jl+1

for each l = 0, · · · , L − 1., i = 1, · · · , n., jl = 0, · · · , dl ., and jl+1 = 1, · · · , dl+1 ..

For simplicity, we introduce a notation r (l) ∈ Rn×dl . such that

(l)
(l)
dxijl (l)
.r
ijl := (l−1)
= 1 − [xijl ]2 , (7.14)
duijl

for all l = 0, · · · , L − 1., i = 1, · · · , n., and jl = 1, · · · , dl .. We can calculate

dl+1 (l) dl+1

(l) ∂J ∂uijl+1
.w
ijl = (l) (l)
= vij(l)l+1 pj(l)
l jl+1
, (7.15)
jl+1 =1 ∂uijl+1 ∂zijl jl+1 =1

for all l = 0, · · · , L − 1., i = 1, · · · , n., and jl = 0, · · · , dl .. Now, we are ready to

state the backward propagation:

(L) (L)
(L−1) ∂J ∂xi1 xi1 − yi (L)
vi1
. = (L) (L−1)
= (L)
{1 − [xi1 ]2 } = Yi − yi , . (7.16)
∂xi1 ∂ui1 1 − [xi1 ]2
(l) (l) (l) (l)
(l−1)
n
∂J ∂zi jl dxijl
n
(l)
nδii − 1 − zijl zi jl (l)
vijl = (l)
= w i jl rijl , (7.17)
i =1 ∂zi jl ∂xij(l)l du(l−1)
ijl i =1 nσj(l)
l
72 7 Batch Normalization

Table 7.1 Batch normalization notations

Notation Dimension Remark Relation
x (l) . n × (dl + 1). Layer x (l+1) = (1n , tanh(u(l) )).
(l) n (l)
μ(l) . 1 × dl . Average μjl = 1
n i=1 xijl .
(l) n (l) (l)
σ (l) . 1 × dl . Standard deviation σjl = 1
n i=1 [xijl − μjl ]2 .
(l) (l)
z(l) . n × (dl + 1). Normalization zi0 = 1; zijl
= [xij(l)l − μ(l) (l)
jl ]/σjl , jl > 0.

p (l) . (dl + 1) × dl+1 . Parameter matrix u(l) = z(l) p (l) .

u(l) . n × dl+1 . Linear combination u(l) = z(l) p (l) .
(l) (l)
r (l) . n × dl+1 . Derivative of activation uijl = 1 − [xijl ]2 , jl > 0.
v (l) . n × dl+1 . Gradient v (l) = ∂J /∂u(l) .
w (l) . n × (dl + 1). Gradient w (l) = ∂J /∂z(l) = v (l) [p(l) ]T .
q (l) . (dl + 1) × dl+1 . Gradient q (l) = ∂J /∂p (l) = [z(l) ]T v (l) .

for each i = 1, · · · , n., l = L − 1, · · · , 1., and jl = 1, · · · , dl .. On account of (7.15),

the above formula provides a backward recurrence relation for v (l) .. Finally, we have

n (l) n
(l) ∂J ∂uijl+1 (l) (l)
.q
jl jl+1 = (l) (l)
= vijl+1 zijl , (7.18)
i=1 uijl+1 ∂pjl jl+1 i=1

for each l = 0, · · · , L − 1., jl = 0, · · · , dl ., and jl+1 = 1, · · · , dl+1 ..

The batch normalization notations are illustrated in Table 7.1.

7.3 Python Code

The Python code for the batch normalization is given below:

Batch Normalization Code

import numpy as np

def append(A):
return [Link](([Link](([Link][0],1)),A),axis=1)
def cut(A):
return A[:,1:]
def ave(A):
n=[Link][0]

(continued)
7.4 Discussions 73

return [Link]([Link]((n,n)),A)/n

def batch_normalization(x,y,hidden_layers,al,iter_max):
n=[Link][0]
d=[[Link][1]-1]+hidden_layers+[1]
L=len(d)-1
p,v,z,mu,s=[None]*L,[None]*L,[None]*L,[None]*L,[None]*L
x=[x]+[None]*L
for l in range(L):
p[l]=[Link](0,1,(d[l]+1,d[l+1]))
for iter in range(iter_max):
for l in range(L):
mu[l]=sum(cut(x[l]))/n
s[l]=[Link](sum((cut(x[l])-mu[l])**2)/n)
z[l]=append((cut(x[l])-mu[l])/s[l])
x[l+1]=append([Link]([Link](z[l],p[l])))
Y=x[L][:,1].reshape(n,1)
v[L-1]=Y-y
for l in range(L-1,0,-1):
w=[Link](v[l],p[l].T)
r=cut(1-x[l]*x[l])
v[l-1]=cut(w-ave(w)-z[l]*ave(z[l]*w))*r/s[l]
for l in range(L):
p[l]=p[l]-al*[Link](z[l].T,v[l])
return Y,p,mu,s

x=[Link]([[1,1,1],[1,1,-1],[1,-1,1],[1,-1,-1]])
y=[Link]([[1],[-1],[-1],[1]])
hidden_layers=[2]
Y,p,mu,s=batch_normalization(x,y,hidden_layers,0.1,1000)
print([Link]((Y,y),axis=1))

7.4 Discussions

In the Python code, we have introduced three auxiliary functions to simplify the
expressions in the main function:
1. The append function is defined as appending a column of ones to the left of
a given matrix. For instance, given an input matrix A ∈ Rn×d ., the output of
append function is the matrix (1n , A) ∈ Rn×(d+1) ..
74 7 Batch Normalization

2. The cut function is defined as deleting the most left column of a given matrix.
It is clear that cut (append(A)) = A., but the equation append(cut (A)) = A. is
not valid unless the first column of A is a column of ones.
3. The ave function is defined as a matrix of identical rows with each row being the
column average of a given matrix. For instance, given A ∈ Rn×d ., the output of
ave function is the matrix (1n , · · · , 1n )A ∈ Rn×d ..

Problems
n n
7.1 Given μ = 1
n i=1 xi . and σ = 1
n i=1 (xi − μ)2 ., prove that

∂μ 1 ∂σ xi − μ
. = , = .
∂xi n ∂xi nσ
n n
7.2 Given μ = 1
n i=1 xi ., σ = 1
n i=1 (xi − μ)2 ., and zi = (xi − μ)/σ ., prove
that
∂zi nδii − 1 − zi zi
. = ,
∂xi nσ

where

1, i=i,
.δii =
0, i i.

7.3 Verify (7.12).

7.4 Verify (7.13).
7.5 Verify (7.16) and (7.17).
Chapter 8
Support Vector Machine

Abstract Given two finite sets of real numbers, S+ . and S− ., where all numbers in
S+ . are positive and all numbers in S− . are negative, one can always select a threshold
value (e.g., 0) that separates the two sets, placing all negative numbers to its left and
all positive numbers to its right. Among infinitely many such thresholds, there exists
an optimal threshold that maximizes the minimal distance from all points in S± . to
the threshold. This optimal threshold forms the foundation of the support vector
machine (SVM). In this chapter, we extend this idea to higher dimensional spaces,
Rd ., and discuss two equivalent optimization formulations for SVMs.

8.1 Notations

Let x ∈ Rn×(d+1) . be the input feature matrix of size n and dimension d, where
xi0 = 1. for i = 1, · · · , n.. Assume that the observation vector is binary y ∈ {1, −1}n .
such that yi ∈ {1, −1}.. For each Xi = (xi0 , · · · , xid )T ∈ Rd+1 ., we denote X̂i =
(xi1 , · · · , xid )T ∈ Rd .. Let p = (p0 , · · · , pd )T ∈ Rd+1 . be the parameter vector.
√ p̂ = (p1 , · · · , pd ) ∈ R .. For any vector u ∈ R ., the norm is
Similarly, we denote T d d

denoted as u uT u.. Given any number t ∈ R., we denote

t, t ≥ 0,
t+ := max{t, 0} =
. (8.1)
0, t ≤ 0.

Given a parameter vector p, the output vector is Y = xp ∈ Rn . such that

d
Yi = (xp)i =
. xij pj = XiT p = p0 + X̂iT p̂.
j =0

The feasible set of parameter vectors in Rd+1 . is defined as

.D := {p ∈ Rd+1 : yi (xp)i ≥ 1, i = 1, · · · , n}. (8.2)

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 75

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
76 8 Support Vector Machine

A parameter vector p ∈ D . is called feasible in the sense that yi (xp)i ≥ 1. for each
i = 1, · · · , n.. It is easily seen that if p ∈ D ., so is cp for any c ≥ 1.. Throughout this
chapter, we assume that the given data {x, y}. is separable; namely, the feasible set
D is non-empty. We also assume that the data is nontrivial: There exist i and j such
that yi = yj .. The normalized feasible set is defined as

E := {p/ p
. p ∈ D}. (8.3)

In other words, E is the set of normalized feasible parameter vectors such that
p 1. and yi (xp)i > 0. for all i = 1, ·, n..
Any parameter vector p ∈ Rd+1 . with p > 0. defines a hyperplane in Rd .:

H (p) := {u ∈ Rd : p0 + uT p̂ = 0}.
. (8.4)

The distance from u ∈ Rd . to the hyperplane H (p). is given by

|p0 + uT p̂|
d(u, H (p)) := min
. u−v . (8.5)
v∈H (p) p

8.2 Optimization Problems

Given a nontrivial separable data x ∈ Rn×(d+1) . and y ∈ {1, −1}n ., we minimize the
objective function

n d
λ
J (p) :=
. C(Yi , yi ) + pj2 , (8.6)
2
i=1 j =1

with cost function

(1 + yi )(1 − Yi )+ + (1 − yi )(1 + Yi )+
C(Yi , yi ) :=
. = (1 − yi Yi )+ , (8.7)
2

where p ∈ Rd+1 . is the parameter vectors, Y = xp ∈ Rn . is the output vector, and

λ > 0. is the regularization parameter. Taking partial derivatives gives

∂J −yi , yi Yi < 1,
. = (8.8)
∂Yi 0, yi Yi > 1,

and
∂J
. =− yi xij + λ(1 − δj 0 )pj , (8.9)
∂pj
yi Yi <1
8.2 Optimization Problems 77

where δj 0 . is the Kronecker delta symbol that equals 1 when j = 0. and 0 when
j 0.. The gradient descent method is given by

(k+1) (k)
pj
. = [1 − αλ(1 − δj 0 )]pj + α yi xij , (8.10)
yi Yi <1

where α > 0. is the learning rate.

Since the data is separable, namely, the feasible set D in (8.2) is non-empty, we
have J (p) = λ p 2 /2. for all p ∈ D ..
As λ → 0+ ., the optimization problem of minimizing J (p). converges to the
following constrained optimization problem.

Constrained Optimization Problem I

Let x ∈ Rn×(d+1) . and y ∈ Rn . be a given nontrivial separable data. Find
p ∈ D ⊂ Rd+1 . such that the objective function

d
J1 (p)
. p pj2 (8.11)
j =1

is minimized.

Proposition 8.1 Assume that the data {x, y}. is nontrivial and separable. The
objective function J1 (p). in (8.11) possesses at least one minimizer, denoted by p ∗ .,
in D .
Proof Since the data is separable, the feasible set D is non-empty. Choose any
q ∈ D ., and denote

d
r :=
. qj2 .
j =1

Since the data is nontrivial, we have r > 0.. Moreover, there exist i1 . and i2 . such that
yi1 = 1. and yi2 = −1.. For any p ∈ D . with p r ., we have p0 + X̂iT1 p̂ ≥ 1. and
p0 + X̂iT2 p̂ ≤ −1.. Consequently,

p0 ≥ 1 − X̂iT1 p̂ ≥ 1 − r X̂i1 ,
.
78 8 Support Vector Machine

and

p0 ≤ −1 − X̂iT2 p̂ ≤ −1 + r X̂i2 .
.

This implies that the subset

Dr := {p ∈ D
. p r}

is bounded. It is easily seen that Dr . is also closed. Hence, the objective function J1 .
possesses a minimizer, denoted by p∗ ., in Dr ., with J1 (p∗ ) ≤ r .. Since J1 (p) > r .
for all p ∈ D \ Dr ., we have J1 (p∗ ) ≤ J1 (p). for all p ∈ D .. This completes the
proof.
Lemma 8.1 For any p ∈ D ., we have

J1 (p) ≥ [min d(X̂i , H (p/ p ))]−1 .

. (8.12)
i

Proof For any p ∈ D . and i = 1, · · · , n., we have yi (p0 + X̂iT p̂) ≥ 1.. Since the data
is nontrivial, we also have p > 0.. Consequently, |p0 +X̂iT p̂| = yi (p0 +X̂iT p̂) ≥ 1.
and

|p0 + X̂iT p̂| 1

d(X̂i , H (p)) =
. ≥ .
p J1 (p)

Since H (p) = H (p/ p )., the inequality (8.12) follows.

Lemma 8.1 suggests us to consider the following constrained optimization
problem.

Constrained Optimization Problem II

Let x ∈ Rn×(d+1) . and y ∈ Rn . be a given nontrivial separable data. Find
p ∈ E ⊂ Rd+1 . such that the objective function

J2 (p) := min |p0 + X̂iT p̂|

. (8.13)
i

is maximized.

Proposition 8.2 Assume that the data {x, y}. is nontrivial and separable. The
objective function J2 (p). in (8.13) possesses at least one minimizer, denoted by p∗∗ .,
in E .
Proof The conclusion follows from the fact that E is compact and J2 . is continuous.
8.3 Equivalence Theorem 79

8.3 Equivalence Theorem

We will prove that the two constrained optimization problems are equivalent.
Theorem 8.1 The constrained optimization problem I has a unique minimizer p∗ ∈
D ., and the constrained optimization problem II has a unique maximizer p ∗∗ ∈ E ..
Moreover, we have

p∗ 1
p∗∗ =
. , J1 (p∗ ) = . (8.14)
p∗ J2 (p∗∗ )

Proof On account of Proposition 8.1 and Proposition 8.2, there exist p∗ ∈ D . and
p ∗∗ ∈ E . such that J1 (p∗ ) ≤ J1 (p). for all p ∈ D . and J2 (p∗∗ ) ≥ J2 (p). for all
p ∈ E .. For any p ∈ E ., we have yi (p0 + X̂iT p̂) > 0. for all i = 1, · · · , n.. Hence,
we obtain

J2 (p) = min{J2+ (p), J2− (p)},

. (8.15)

where

.J2+ (p) = min (p0 + X̂iT p̂) = p0 + min (X̂iT p̂), . (8.16)
yi =1 yi =1

J2− (p) = min (−p0 − X̂iT p̂) = −p0 + min (−X̂iT p̂). (8.17)
yi =−1 yi =−1

Note that J2 (p∗∗ ) ≥ J2 (p). for all p ∈ E ..

We claim J2+ (p∗∗ ) = J2− (p∗∗ ).; namely, p0∗∗ = g(p̂∗∗ ). with

1
g(u) :=
. [ min (−X̂iT u) − min (X̂iT u)], u ∈ Rd .
2 yi =−1 yi =1

Otherwise, we may choose p ∈ E . with p̂ = p̂ ∗∗ . and p0 = g(p̂). to obtain

J2+ (p) + J2− (p) J + (p∗∗ ) + J2− (p∗∗ )

J2 (p) = J2+ (p) = J2− (p) =
. = 2 > J2 (p∗∗ ),
2 2
a contradiction.
Now, we define a function f ∈ C(Rd , R). as

1
.f (u) := [ min (−X̂iT u) + min (X̂iT u)], (8.18)
2 yi =−1 yi =1

for u ∈ Rd .. It follows from the above argument that J2 (p∗∗ ) = f (p̂ ∗∗ ).. Moreover,
by definition, we observe that
80 8 Support Vector Machine

f (cu) = cf (u), f (u + v) ≥ f (u) + f (v),

. (8.19)

for any c > 0. and u, v ∈ Rd .. Since p 1. and f (p̂) > 0. for p ∈ E ., the function
f (u). has a unique maximizer on the unit sphere S := {u ∈ Rd , u 1}..
We show that the maximizer of f on S is exactly p̂∗∗ .. This is because for any
u ∈ S . with f (u) > 0., we may choose p ∈ Rd+1 . with p̂ = u. and p0 = g(u).
such that p ∈ E . and J2 (p) = J2+ (p) = J2− (p) = f (p̂)., and hence f (p̂ ∗∗ ) =
J2 (p∗∗ ) ≥ J2 (p) = f (p̂) = f (u).. On the other hand, if u ∈ S . with f (u) ≤ 0., then
f (p̂ ∗∗ ) > 0 ≥ f (u).. This proves that p̂ ∗∗ . is the unique maximizer of f on S. In
view of p0∗∗ = g(p̂∗∗ )., the maximizer of J2 . is also unique.
Finally, since p ∗ /J1 (p∗ ) ∈ E . and p∗∗ /J2 (p∗∗ ) ∈ D ., we have

1 1
J2 (p∗∗ ) ≥ J2 (p∗ /J1 (p∗ )) ≥
. ≥ = J2 (p∗∗ ).
J1 (p∗ ) J1 (p∗∗ /J2 (p∗∗ ))

The above inequalities should be equalities, and the identities (8.14) are proved.
Finally, the uniqueness of p ∗ . follows from the uniqueness of p∗∗ ..

Problems
8.1 Given p ∈ Rd+1 . with p > 0.. Let H (p). be the hyperplane defined in (8.4).
For any u ∈ Rd ., prove that

|p0 + uTp̂|
. min u−v .
v∈H (p) p
8.2 Verify that H (p) = H (cp). for any c 0. and p ∈ Rd+1 . with p > 0..
8.3 Given Yi ∈ R. and yi ∈ {1, −1}., prove that

(1 + yi )(1 − Yi )+ + (1 − yi )(1 + Yi )+
. = (1 − yi Yi )+ .
2

8.4 Assume that D in (8.2) is non-empty. Prove that J (p) = λ p 2 /2. for all p ∈
D ..
8.5 Prove that the data {x, y}. is nontrivial if and only if p > 0. for any p ∈ D ..
8.6 Assume that the data {x, y}. is nontrivial and separable. Prove that the feasible
set D defined in (8.2) is closed in Rd+1 ..
8.7 Assume that the data {x, y}. is nontrivial and separable. Prove that the normal-
ized feasible set E defined in (8.3) is compact (i.e., bounded and closed) in Rd+1 ..
8.8 Show that p/J2 (p) ∈ D . if p ∈ E ..
8.9 Prove that E defined in (8.3) is the same as

E = {p ∈ Rd+1
. p 1, yi (xp)i > 0, i = 1, · · · , n}.
8.3 Equivalence Theorem 81

8.10 Let E and f be defined as in (8.3) and (8.18), respectively. Prove that f (p̂) >
0. for any p ∈ E ..
8.11 Verify (8.19).
8.12 Assume h ∈ C(Rd , R). is convex and h(cu) = ch(u). for any c > 0. and
u ∈ Rd .. If the minimum of h on the unit sphere S = {u ∈ Rd u 1}. is
negative, prove that the minimizer on S is unique.
Chapter 9
Gradient Methods

Abstract This chapter introduces several fundamental gradient-based iterative

methods for solving linear systems and optimization problems, including the
conjugate gradient method, the Jacobi method, and the heavy ball method. We
begin with the conjugate gradient method for symmetric positive definite matrices,
presenting the derivation of the iteration formulas, orthogonality properties, and the
underlying vector space structure. The method is shown to converge to the exact
solution in at most m steps for an m-dimensional system. Next, we discuss one-step
iterative methods, including the Jacobi and gradient descent methods, and analyze
their convergence using eigenvalue techniques. The concept of a learning rate is
introduced, and we formulate a minimax optimization problem to determine the
optimal learning rate for the fastest convergence. Finally, the chapter extends these
ideas to two-step methods, or heavy ball methods, which incorporate a momentum
term to accelerate convergence. We analyze the spectral radius of the associated
iteration matrix and present a minimax problem to optimize both the learning rate
and momentum parameter.

9.1 Conjugate Gradient Method

Let A ∈ Rm×m . be a symmetric and positive definite matrix. We may define an inner
product on Rm . as

(u, w) = uT Aw, u, w ∈ Rm .
. (9.1)

Given b ∈ Rm ., we intend to solve the linear system

Av = b.
. (9.2)

The above equation is a generalization of (2.12) in linear regression. In general, the

computation cost of v = A−1 b. is O(m3 ).. We shall introduce an iteration method to
find an approximate solution using less computation time.

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 83

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
84 9 Gradient Methods

Let v1 ∈ Rm . be any initial guess of solution, and we denote e1 = v1 − v . and

r1 = Av1 − b = Ae1 . to be the error vector and residual vector, respectively. We also
initialize u1 = r1 .. If u1 = 0., then we have found the solution. Otherwise, we define
r2 = Ae2 . with

(e1 , u1 )
e2 = e1 − α1 u1 , α1 =
. .
(u1 , u1 )

Note that α1 u1 . is the projection of e1 . on u1 . and e2 . is perpendicular to u1 .; namely,

(e2 , u1 ) = 0.. We also define

(r2 , u1 )
u2 = r2 −
. u1
(u1 , u1 )

such that (u2 , u1 ) = 0..

Now, we assume that ej ., rj = Aej . and uj . with j = 1, · · · , k . have been
calculated. If uk = 0., then we are done. Otherwise, we define

rk+1 = Aek+1 ,
. (9.3)

with
(ek , uk )
ek+1 = ek − αk uk , αk =
. . (9.4)
(uk , uk )

We also use the Gram-Schmidt process to define

k
(rk , uj )
.uk+1 = rk+1 − uj (9.5)
(uj , uj )
j =1

such that

. (uk , uj ) = 0, j < k. (9.6)

We have the following proposition.

Proposition 9.1 The following four vector spaces are identical:

. Sk(1) := span{r1 , r2 , · · · , rk },
(2)
Sk := span{u1 , u2 , · · · , uk },
(3)
Sk := span{r1 , Ar1 , · · · , Ark−1 },
(4)
Sk := span{u1 , Au1 , · · · , Auk−1 }.
9.1 Conjugate Gradient Method 85

Proof We prove the statement by induction. First, since u1 = r1 ., we have

(1) (2) (3) (4)

Sk = Sk = Sk = Sk
.

for k = 1.. Now, we assume that the above identities are true for all k = 1, · · · , K ..
(1) (2) (3) (4)
In view of (9.5), we have Sk = Sk . and Sk = Sk . for k = K + 1.. Moreover,
since

rk+1 − rk = A(ek+1 − ek ) = −αk Auk

.

by (9.4), we also have Sk(1) = Sk(4) . for k = K + 1.. Therefore, the four vector spaces
are identical for k = K + 1.. This completes the proof.
On account of Proposition 9.1, we denote

Sk := span{r1 , r2 , · · · , rk } = span{u1 , u2 , · · · , uk }
.

= span{r1 , Ar1 , · · · , Ark−1 } = span{u1 , Au1 , · · · , Auk−1 }. (9.7)

We have the following corollary.

Corollary 9.1 For any j < k ., we have

rkT uj = (ek , uj ) = 0, .
. (9.8)
(rk+1 , uj ) = (ek+1 , Auj ) = 0. (9.9)

Proof If k = 2., then (9.4) implies that r2T u1 = (e2 , u1 ) = 0.. Now, we assume
that (9.8) holds for any 1 ≤ j < k ≤ K .. By (9.6) and (9.4), we have (ek+1 , uk ) = 0.
and

(ek+1 , uj ) = (ek , uj ) − αk (uk , uj ) = 0

.

for 1 ≤ j < k = K .. Thus, (9.8) is valid for all 1 ≤ j < k = K + 1.. By

induction, (9.8) is proved for all 1 ≤ j < k .. Hence, ek+1 . is perpendicular to the
space Sk . defined in (9.7), which yields (9.9).
In view of (9.9), we can rewrite (9.5) as

uk+1 = rk+1 + βk uk ,
. (9.10)

where

(rk+1 , uk ) r T (rk+1 − rk ) T r
rk+1 k+1
.βk := − = k+1 = . (9.11)
(uk , uk ) uk (rk − rk+1 )
T T
rk r k

The projection (9.4) can be rewritten as

86 9 Gradient Methods

vk+1 = vk − αk uk ,
. (9.12)

where

(ek , uk ) r T rk
αk =
. = Tk . (9.13)
(uk , uk ) uk Auk

Conjugate Gradient Method Pseudo Code

Input: A ∈ Rm×m , b ∈ Rm .
Initialize: v1 = 0 ∈ Rm , r1 = Av1 − b, u1 = r1 .
Loop: for k = 1, 2, · · · ., while rkT rk > ε., compute the following two lines
• Projection: vk+1 = vk − αk uk ., rk+1 = rk − αk Auk . where αk =
rkT rk /(uTk Auk ).
• Orthogonalization: uk+1 = rk+1 + βk uk . where βk = rk+1
T r T
k+1 /(rk rk ).

Output: vk .

Theorem 9.1 If u1 , · · · , um . are nonzero vectors, then rm+1 = 0., em+1 = 0., and
vm+1 = v ..
Proof Since (uk , uj ) = 0. for j < k ., the nonzero vectors u1 , · · · , um . are
orthogonal to each other. Hence, Sm = span{u1 , · · · , um } = Rm . and u1 , · · · , um .
are orthogonal basis in Rm .. By (9.9), we have (rm+1 , uj ) = 0. for all j = 1, · · · , m..
Consequently, rm+1 = 0., which implies em+1 = A−1 rm+1 = 0. and vm+1 =
v + em+1 = v .. This completes the proof.

9.2 One-Step Method

Let A ∈ Rm×m . be a nonsingular matrix. Choose a diagonal matrix D ∈ Rm×m .

with positive diagonal entries such that all eigenvalues of D −1 A. are positive. The
linear system Av = b. can be rewritten as Dv = Dv + α(b − Av). for any α > 0..
Multiplying by D −1 . on both sides gives

v = D −1 [Dv + α(b − Av)] = v − αD −1 (Av − b).

.

This motivates us to consider the one-step iterative method

.vk+1 = vk − αD −1 (Avk − b). (9.14)

9.2 One-Step Method 87

Here, α > 0. is the learning rate. If α = 1. and Djj = Ajj . for j = 1, · · · , m., then
the above iteration is called the Jacobi method. If D = I . and A is a positive definite
symmetric matrix, then the above iteration is the gradient descent method for the
optimization problem of minimizing the quadratic function

1 T
Q(v) :=
. v Av − v T b.
2
Let ek = vk − v . be the error vector. The iteration (9.14) is rewritten as

ek+1 = (I − αD −1 A)ek = (I − αD −1 A)k e1 .

. (9.15)

Let 0 < λ1 ≤ · · · ≤ λm . be the eigenvalues of D −1 A. (counting multiplicity). Denote

J1 (α) := ρ(I − αD −1 A) =
. max |1 − αλj |. (9.16)
j =1,··· ,m

If α > 0. is sufficiently small, say α < 1/λm ., then J1 (α) < 1. and ek → 0. as k →
∞.. To obtain a fast convergence rate, we shall investigate the following minimax
optimization problem.

Minimax Optimization Problem I

Let 0 < λ1 ≤ · · · ≤ λm . be the eigenvalues of D −1 A. (counting multiplicity).
Find α > 0. such that the objective function J1 (α). defined in (9.16) is
minimized.

Theorem 9.2 Let α ∗ = 2/(λm + λ1 ).. We have

λ m − λ1
J1 (α) > J1 (α ∗ ) =
. (9.17)
λ m + λ1

for all α > 0. with α α ∗ ..

Proof Given α > 0., then function fα (λ) = |1 − αλ|. is a convex function, and
its maximum value on [λ1 , λm ]. is achieved at either λ1 . or λm .. Hence, we obtain
from (9.16) that

J1 (α) = max{fα (λ1 ), fα (λm )}.

.

If α = α ∗ = 2/(λm + λ1 )., then fα (λ1 ) = fα (λm ) = (λm − λ1 )/(λm + λ1 ).. This

implies that

J1 (α ∗ ) = (λm − λ1 )/(λm + λ1 ).
.
88 9 Gradient Methods

If α > α ∗ ., then

J1 (α) ≥ fα (λm ) = αλm − 1 > α ∗ λm − 1 = J1 (α ∗ ).

.

If α ∈ (0, α ∗ )., then

J1 (α) ≥ fα (λ1 ) = 1 − αλ1 − 1 > 1 − α ∗ λ1 = J1 (α ∗ ).

.

This completes the proof.

9.3 Two-Step Method

We generalize the one-step method (9.14) as

vk+1 = vk − αD −1 (Avk − b) + β(vk − vk−1 ),

. (9.18)

where α > 0. is the learning rate, and β ≥ 0. is another hyperparameter related to

the moment. When β = 0., the above two-step iterative method reduces to the one-
step method in (9.14). In the literature, this two-step method is also called the heavy
ball method. Now, the error vector ek = vk − v . satisfies a second-order difference
equation:

ek+1 = [(1 + β)I − αD −1 A]ek − βek−1 ,

. (9.19)

which can be rewritten as

ek+1 (1 + β)I − αD −1 A − βI ek
. = . (9.20)
ek I 0 ek−1

Let 0 < λ1 ≤ · · · ≤ λm . be the eigenvalues of D −1 A.; then the eigenvalues of the

block matrix

(1 + β)I − αD −1 A − βI
.
I 0

are

1 + β − αλj ± (1 + β − αλj )2 − 4β
. , j = 1, · · · , m.
2
This motivates us to consider the following minimax optimization problem.
9.3 Two-Step Method 89

Minimax Optimization Problem II

Let 0 < λ1 ≤ · · · ≤ λm . be the eigenvalues of D −1 A. (counting multiplicity).
Find α > 0. and β ≥ 0. such that the objective function

J2 (α, β) :=
. max ρj (α, β) (9.21)
j =1,··· ,m

is minimized, where

|1 + β − αλj | + (1 + β − αλj )2 − 4β
ρj (α, β) :=
. (9.22)
2

is the spectral radius of the matrix

(1 + β)I − αλj −β
. .
1 0

Theorem 9.3 Let

2 √ √ 2
∗ 2 ∗ λm − λ1
.α = √ √ , β = √ √ . (9.23)
λm + λ1 λm + λ1

We have
√ √
λm − λ1
J2 (α, β) > J2 (α ∗ , β ∗ ) = √
. √ (9.24)
λm + λ1

for all α > 0. and β ≥ 0. with (α, β) (α ∗ , β ∗ )..

Proof It is easy to verify that
√ √
1 + β ∗ − α ∗ λ1 1 + β ∗ − α ∗ λm λm − λ1
. =− = β∗ =√ √ ,
2 2 λm + λ1

and
√ √
λm − λ1
.J2 (α ∗ , β ∗ ) = ρ1 (α ∗ , β ∗ ) = ρm (α ∗ , β ∗ ) = √ √ .
λm + λ1

We consider the following four cases, respectively:

(i) If 0 < α < α ∗ . and β ≥ β ∗ ., then we have
90 9 Gradient Methods

1 + β − αλ1 1 + β ∗ − α ∗ λ1
J2 (α, β) ≥ρ1 (α, β) ≥
. > = J2 (α ∗ , β ∗ ).
2 2

(ii) If α > α ∗ . and 0 ≤ β ≤ β ∗ ., then we have

1 + β − αλm 1 + β ∗ − α ∗ λm
J2 (α, β) ≥ρm (α, β) ≥ −
. >− = J2 (α ∗ , β ∗ ).
2 2

(iii) If 0 < α ≤ α ∗ . and 0 ≤ β < β ∗ ., then we have (1 + β − αλ1 )2 − 4β > 0. and

1 + β − αλ1 + (1 + β − αλ1 )2 − 4β
J2 (α, β) ≥ρ1 (α, β) =
.
2
1 + β − α ∗ λ1 + (1 + β − α ∗ λ1 )2 − 4β
≥
2
1 + β ∗ − α ∗ λ1 + (1 + β ∗ − α ∗ λ1 )2 − 4β ∗
> = J2 (α ∗ , β ∗ ).
2

(iv) If α ≥ α ∗ . and β > β ∗ ., then we have

.J2 (α, β) ≥ρm (α, β) ≥ β> β ∗ = J2 (α ∗ , β ∗ ).

A combination of the above four cases yields the desired results.

Problems
9.1 Let rk ., ek ., and uk . be defined recursively from (9.3), (9.4), and (9.5). Prove that
T r = 0., r T u = 0., and uT r = r T r ..
rk+1 k k+1 k k k k k
9.2 Verify (9.11) and (9.13).
9.3 Prove that if the iteration (9.14) converges, then the limit is a solution to the
linear system Av = b..
9.4 Verify (9.15).
9.5 Given α > 0. and 0 < λ1 < λm ., prove that the function fα (λ) = |1 − αλ|. is a
convex function, and its maximum value on [λ1 , λm ]. is achieved at either λ1 . or λm ..
9.6 Given x ∈ Rn×(d+1) . and y ∈ Rn . with det(x T x) 0., show that the minimizer
of the error SSE xp − y 22 . can be approximated by the iteration (9.14) where
D = I ., A = x T x ., b = x T y ., and the learning rate α > 0. is chosen to be sufficiently
small.
9.7 Let 0 < λ1 ≤ · · · ≤ λm . be the eigenvalues of D −1 A.; prove that the eigenvalues
of the block matrix
9.3 Two-Step Method 91

(1 + β)I − αD −1 A − βI
.
I 0

are

1 + β − αλj ± (1 + β − αλj )2 − 4β
. , k = 1, · · · , m.
2
9.8 Given α > 0., β ≥ 0., and λj > 0., show that ρj (α, β). defined in (9.22) is the
spectral radius of the matrix

(1 + β)I − αλj −β
. .
1 0

9.9 Given 0 < λ1 < λm ., let α ∗ . and β ∗ . be defined as in (9.23). For any λj ∈
(λ1 , λm )., prove that
√ √
λm − λ1
ρj (α ∗ , β ∗ ) < ρ1 (α ∗ , β ∗ ) = ρm (α ∗ , β ∗ ) = √
. √ ,
λm + λ1

where ρj (α, β). is defined as in (9.22).

9.10 Given 0 < λ1 < λm ., let α ∗ . and β ∗ . be defined as in (9.23). Prove that the
function

f (β) := 1 + β − α ∗ λ1 +
. (1 + β − α ∗ λ1 )2 − 4β

is strictly decreasing (i.e., f (β) < 0.) for β ∈ [0, β ∗ ]..

√
9.11 For any β ≥ 0. and γ ≥ 0., prove that |γ + γ 2 − 4β| ≥ 2 β ..
Chapter 10
Dimensionality Reduction

Abstract This chapter presents the theory of singular value decomposition (SVD)
and principal component analysis (PCA), which are fundamental tools for ana-
lyzing and approximating high-dimensional data. We begin by reviewing the
Schur decomposition for square matrices and unitary transformations, laying the
groundwork for understanding SVD. The chapter then introduces singular value
decomposition, proving its existence, uniqueness, and the relationship between
singular values and the eigenvalues of AT A. and AAT .. Building on SVD, we develop
principal component analysis as a method for dimensionality reduction. We define
the Rayleigh quotient and formulate a sequence of optimization problems whose
solutions yield the principal components of a given matrix. The chapter establishes
the orthogonality of left and right singular vectors and demonstrates that the sum
of the first k principal components provides the best rank-k approximation of the
original matrix in the Frobenius norm.

10.1 Schur Decomposition

Recall that a matrix U ∈ Cm×m . is unitary if U U H = U H U = I ., where U H ∈

Cm×m . is the conjugate transpose of U defined as UjHk = Ukj ..
Theorem 10.1 For any A ∈ Cm×m ., there exist a unitary matrix U ∈ Cm×m . (i.e.,
U U H = U H U = I .) and an upper triangular matrix T ∈ Cm×m . (i.e., Tj k = 0. for
1 ≤ k < j ≤ m.) such that A = U T U H .. We may order the diagonal terms of T
such that |T11 | ≥ |T22 | ≥ · · · ≥ |Tmm |..
Proof We prove by induction. The result is trivial when m = 1.. For a general m >
1., assume the result is true for m − 1.. Let λ1 . be the eigenvalue of A with the largest
modulus. Let u1 ∈ Cm . be the associated eigenvector with u1 2 = 1.. By a standard
Gram-Schmidt process, we can find m − 1. orthonormal vectors u2 , · · · , um ∈ Cm .
such that U1 = (u1 , · · · , um ) ∈ Cm×m . is unitary. Since Au1 = λ1 u1 ., we have

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 93

X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
94 10 Dimensionality Reduction

λ1 α1H
AU1 = U1
. ,
0 A2

where α1 ∈ Cm−1 . and A2 ∈ C(m−1)×(m−1) .. By inductive assumption, there

exist a unitary matrix U2 ∈ C(m−1)×(m−1) . and an upper triangular matrix T2 ∈
C(m−1)×(m−1) . such that A2 = U2 T2 U2H .. Substituting this into the above equation
gives

λ1 α1H
A = U1
. U1H .
0 U2 T2 U2H

In view of the matrix factorization

λ1 α1H 1 0 λ1 α1H U2 1 0
. = ,
0 U2 T2 U2H 0 U2 0 T2 0 U2H

we obtain A = U T U H ., where

1 0
U=
.
0 U2

is unitary and

λ1 α1H U2
T =
.
0 T2

is upper triangular. This completes the proof.

10.2 Singular Value Decomposition

Theorem 10.2 For any matrix A ∈ Cm×n ., there exist a unitary matrix U ∈ Cm×m .,
a unitary matrix V ∈ Cn×n ., and a rectangular diagonal matrix Σ ∈ Cm×n . with
nonnegative diagonal terms, such that A = U ΣV H .. We may order the diagonal
terms of Σ . such that Σjj ≥ Σkk . for 1 ≤ j < k ≤ min{m, n}..
Proof Without loss of generality, we may assume that m ≤ n. and A is of rank m;
namely, AAH ∈ Cm×m . is a positive definite Hermitian matrix. There exist a unitary
matrix U ∈ Cm×m . and a diagonal matrix D ∈ Cm×m . with positive diagonal terms
such that AAH = U D 2 U H .. In particular, we may order the diagonal terms of D as
D11 ≥ D22 ≥ · · · ≥ Dmm > 0.. Set V1 = AH U D −1 ∈ Cn×m .. It is easily seen that
A = U DV1H . and
10.3 Principal Component Analysis 95

.V1H V1 = D −1 U H AAH U H D −1 = D −1 D 2 D −1 = I.

Finally, we set Σ = (D, 0) ∈ Rm×n . and use a Gram-Schmidt process to choose

V2 ∈ Cn×(n−m) . such that V H = (V1H , V2H ) ∈ Cn×n . is unitary. This completes the
proof.
The decomposition A = U ΣV H . is called a singular value decomposition, and
the nonnegative diagonal terms of Σ . are called the singular values of A.
Remark 10.1 If A = U ΣV H ., then AAH = U ΣΣ H U H . and AH A = V Σ H ΣV H ..
Therefore, the square of any singular value of A is an eigenvalue of AAH . (or AH A.).
In particular, if A is Hermitian, then the absolute values of the eigenvalues of A are
the singular values of A. However, for a general square matrix A, the absolute value
of eigenvalues and the singular values may be different.
The following proposition shows that the rectangular diagonal matrix of the
singular value decomposition is unique up to a permutation of the singular values.
Proposition 10.1 Assume A = U ΣV H = Ũ Σ̃ Ṽ H . are two singular value
decompositions of A ∈ Cm×n ., where U and Ũ . are unitary matrices in Cm×m ., V and
Ṽ . are unitary matrices in Cn×n ., and Σ . and Σ̃ . are rectangular diagonal matrices in
Cm×n .. If the diagonal terms of Σ . and Σ̃ . are ordered as Σjj ≥ Σkk . and Σ̃jj ≥ Σ̃kk .
for 1 ≤ j < k ≤ min{m, n}., then we have Σ = Σ̃ ..
Proof Without loss of generality, we may assume m ≤ n. and rewrite Σ =
(D, 0). and Σ̃ = (D̃, 0)., where D and D̃ . are diagonal matrices in Cm×m .. Since
AAH = U D 2 U H = Ũ D̃ 2 Ũ H ., the diagonal terms of D 2 . and D̃ 2 . are the same
as the eigenvalues of AAH .. By assumption, these diagonal terms are ordered as
D11 ≥ · · · ≥ Dmm . and D̃11 ≥ · · · ≥ D̃mm .. Therefore, we have D = D̃ . and
Σ = Σ̃ .. This completes the proof.

10.3 Principal Component Analysis

Throughout this section, we assume that a matrix A ∈ Rm×n . is given. By Schur

decomposition, we obtain AT A = V ΛV T ., where V = (v1 , · · · , vn ). is an
orthogonal matrix (i.e., V V T = V T V = I .) and Λ. is a diagonal matrix with
nonnegative diagonal terms λ1 ≥ λ2 ≥ · · · ≥ λn ≥ 0..
The first Rayleigh quotient is defined as

v T AT Av Av 22
R1 (v) :=
. = , (10.1)
vT v v 22

for v ∈ Rn . with v 0..

96 10 Dimensionality Reduction

The First Optimization Problem

Find v ∈ Rn . with v 2 = 1. such that R1 (v) Av 2 . is
2 maximized.

Proposition 10.2 For any normalized vector v ∈ Rn ., we have R1 (v) ≤ R1 (v1 ) =

λ1 ..
Proof For any v = c1 v1 + · · · + cn vn ∈ Rn . with c12 + · · · + cn2 = 1., we have

n n
.R1 (v) = v T AT Av = λj cj2 ≤ λ1 cj2 = λ1 = R1 (v1 ).
j =1 j =1

This completes the proof.

√
For each k = 1, · · · , min{m, n}., we denote σk = λk . to be the k-th singular value.
The right eigenvector vk . of AT A. associated with the eigenvalue λk . is called the right
singular vector of A associated with σk .. If σk > 0., we also define uk = (Avk )/σk ∈
Rm . as the left singular vector of A associated with σk ..
Lemma 10.1 For each k = 1, · · · , min{m, n}. with σk > 0., the left singular vectors
u1 , · · · , uk . are orthonormal vectors in Rm ..
Proof For any j = 1, · · · , k ., by definition, we have uTj uj = (vjT AT Avj )/λj = 1.,
and hence uj . is a normalized vector. If 1 ≤ j < k ≤ min{m, n}. and σk > 0.,
then uTj uk = (vjT AT Avk )/(σj σk ) = 0.. Therefore, by induction, for each k =
1, · · · , min{m, n}. with σk > 0., the vectors u1 , · · · , uk . are orthonormal vectors.
For each k = 1, · · · , min{m, n}., we define the k-th principal component of A as

Pk := Avk vkT .
. (10.2)

We can regard Pk . as the projection of A on vk .. The k-th Rayleigh quotient is given

by

v T ATk Ak v Ak v 22
Rk (v) :=
. = , (10.3)
T
v v v 22

for any nonzero vector v ∈ Rn ., where

k−1 k−1
Ak := A −
. Pj = A − Avj vjT . (10.4)
j =1 j =1

If v ∈ Rn . satisfies v 2 = 1. and v T vj = 0. for j = 1, · · · , k − 1., then Ak v = Av .

and Rk (v) Ak v 22 Av 22 ..
10.3 Principal Component Analysis 97

The k-th Optimization Problem

For each k = 2, · · · , min{m, n}., find v ∈ Rn . with v 2 = 1. and v T vj = 0. for
j = 1, · · · , k − 1., such that Rk (v) Ak v 22 Av 22 . defined as in (10.3)
is maximized.

Proposition 10.3 For each k = 2, · · · , min{m, n}., we have Rk (v) ≤ Rk (vk ) = λk .

for any normalized vector v ∈ Rn . with v T vj = 0. for j = 1, · · · , k − 1..
Proof For any v = c1 v1 + · · · + cn vn ∈ Rn . with c12 + · · · + cn2 = 1. and v T vj = 0.
for j = 1, · · · , k − 1., we have c1 = · · · = ck−1 = 0., and
n n
Rk (v) = v T AT Av =
. λj cj2 ≤ λk cj2 = λk = Rk (vk ).
j =k j =k

This completes the proof.

Theorem 10.3 If σk = 0. for some k = 1, · · · , min{m, n}., then Ak = 0.; namely,

k−1
A=
. Pj . (10.5)
j =1

If σk > 0. for all k = 1, · · · , min{m, n}., then we have

min{m,n}
A=
. Pj . (10.6)
j =1

Proof From (10.4), we have Ak vj = 0. for any 1 ≤ j < k .. If σk = 0. for some k =

1, · · · , min{m, n}., then Proposition 10.3 implies that Ak vj = 0. for all k ≤ j ≤ n..
Therefore, Ak v = 0. for any v = c1 v2 + · · · + cn vn ∈ Rn .. Consequently, Ak = 0.,
and (10.5) is proved.
Now, we assume that σk > 0. for all k = 1, · · · , min{m, n}.. If m ≥ n., then
by (10.4), we obtain An+1 vj = 0. for all j = 1, · · · , n., which implies An+1 = 0.. If
m ≤ n., then from uj = (Auj )/σj . we obtain σj uTj Am+1 = vjT AT Am+1 = 0. for all
j = 1, · · · , m., which implies that Am+1 = 0.. This proves (10.6).
The following lemma is a generalization of triangle inequality and Pythagorean
theorem.
Lemma 10.2 Let w1 , · · · , wk ∈ Rn . be an orthonormal basis of a subspace W ⊂
Rn .. For any u ∈ Rn . and w ∈ W ., we have
98 10 Dimensionality Reduction

k
. u−w 2
2 u 2
2 − (uT wj )2 ,
j =1

where the equality holds if and only if

k
.w= (uT wj )wj
j =1

is the projection of u on W .
Proof If k = 1., then

. u − cw1 2
2 = c2 − 2(uT w1 )c + uT u u 2
2 − (uT w1 )2 ,

where the equality holds if and only if c = uT w1 .. If k ≥ 2., by induction, we obtain

k k−1 k
. u− cj wj 2
2 u− cj wj 2
2 − (uT wk )2 u 2
2 − (uT wj )2 ,
j =1 j =1 j =1

where the equality holds if and only if cj = uT wj . for all j = 1, · · · , k .. This

completes the proof.
Lemma 10.3 Let W be a subspace of Rn . with dimension k ≤ n.. For any vectors
e1 , · · · , ek−1 ∈ Rn ., there exists w ∈ W . such that w 2 = 1. and w T ej = 0. for all
j = 1, · · · , k − 1..
Proof Let w1 , · · · , wk . be an orthonormal basis of W . We define a matrix

B := (e1 , · · · , ek−1 )T (w1 , · · · , wk ) ∈ R(k−1)×k ,

.

namely, Bj l = ejT wl . for j = 1, · · · , k − 1. and l = 1, · · · , k .. The column vectors of

. B T = (w1 , · · · , wk )T (e1 , · · · , ek−1 ) ∈ R(k−1)×k

are denoted by b1 , · · · , bk−1 ∈ Rk .; namely, B T = (b1 , · · · , bk ). with

bj = (w1 , · · · , wk )T ej , j = 1, · · · , k − 1.
.

Since the subspace spanned by these vectors has a dimension less than k, there
exists a vector c = (c1 , · · · , ck )T ∈ Rk . such that c 2 = 1. and cT bj = 0. for all
j = 1, · · · , k − 1.. This implies that Bc = 0.. By choosing

.w = w1 c1 + · · · + wk ck ∈ Rn ,
10.3 Principal Component Analysis 99

we obtain

ejT w = Bj 1 c1 + · · · + Bj k ck = (Bc)j = 0
.

for each j = 1, · · · , k − 1.. Moreover, w 2

2 = c12 + · · · + ck2 = 1.. The proof is
completed.
The following theorem demonstrates that the sum of first k principal components
provides the best approximation to A among all matrices of rank k with respect to
the Frobenius norm. Recall that the Frobenius norm of B ∈ Rm×n . is the same as the
l 2 .-norm of B when viewed as a vector in Rmn .; namely,

m n
. B F := |Bj k |2 . (10.7)
j =1 k=1

Theorem 10.4 For any B ∈ Rm×n . with rank(B) = k ., we have

n
. A−B 2
F ≥ λj ,
j =k+1

where the equality holds if B = P1 + · · · + Pk . is the sum of first k principal

components.
Proof Denote the column vectors of AT . and B T . by a1 , · · · , am ∈ Rn . and
b1 , · · · , bm ∈ Rn ., respectively. Then we have
m
. A−B 2
F = al − bl 2
2.
l=1

Since rank(B) = k ., the subspace W = span{b1 , · · · , bm }. has a dimension k. Recall

that v1 , · · · , vn . are eigenvectors of AT A. corresponding to λ1 , · · · , λn ., respectively.
These vectors form an orthonormal basis in Rn .. By repeated applications of
Lemma 10.3, we can find an orthonormal basis w1 , · · · , wk . in W such that wjT vl =
0. for any 1 ≤ l < j ≤ k .. For each l = 1, · · · , m., it follows from Lemma 10.2 that

k
. al − bl 2
2 al 2
2 − (alT wj )2 ,
j =1

where the equality holds if and only if

k
bl =
. (alT wj )wj
j =1
100 10 Dimensionality Reduction

is the projection of al . on W . Hence, we obtain

m m k k
. A−B 2
F ≥ al 2
2 − (alT wj )2 A 2
F − Awj 2
2.
l=1 l=1 j =1 j =1

For each j = 1, · · · , k ., on account of Propositions 10.2 and 10.3, we have

Awj 22 Avj 22 = λj .. Consequently, we obtain

n k n
. A−B 2
F ≥ λj − λj = λj ,
j =1 j =1 j =k+1

where the equality holds if

k
bl =
. (alT vj )vj ;
j =1

namely,

k k
B=
. Avj vjT = Pj .
j =1 j =1

This completes the proof.

Problems
10.1 For any u1 ∈ Cn . with u1 2 = 1., show that there exist normalized
orthonormal vectors u2 , · · · , um ∈ Cm . such that U1 = (u1 , · · · , um ) ∈ Cm×m .
is unitary.
10.2 Verify

λ1 α1H 1 0 λ1 α1H U2 1 0
. = .
0 U2 T2 U2H 0 U2 0 T2 0 U2H

10.3 If A ∈ Cm×m . is normal, i.e., AAH = AH A., prove that A is unitarily

diagonalizable; namely, A = U DU H . for some unitary matrix U and diagonal
matrix D.
10.4 If A ∈ Rm×m . is real symmetric (i.e., AT = A.), prove that there exist an
orthogonal matrix V ∈ Rm×m . (i.e., V V T = V T V = I .) and a diagonal matrix
D ∈ Rm×m . such that A = V DV T ..
10.3 Principal Component Analysis 101

10.5 Assume m ≤ n.. For any V1 ∈ Cn×m . with V1H V1 = I ., prove that there exists
V2 ∈ Cn×(n−m) . such that V H = (V1H , V2H ) ∈ Cn×n . is unitary.
10.6 Compute the eigenvalues and singular values of the matrix

12
A=
. .
01

10.7 Let A ∈ Rm×n . with the first Rayleigh quotient R1 (v). defined as in (10.1).
For any nonzero vector v ∈ Rn . and its normalization w = v/ v 2 .. Prove that
R1 (v) = R1 (w) Aw 22 ..
10.8 Given A ∈ Rm×n ., let Ak ∈ Rm×n . defined as in (10.4). Prove that

k−1
ATk Ak = AT A −
. λj vj vjT .
j =1

10.9 Given any vectors b1 , · · · , bk−1 ∈ Rk ., prove that there exists a vector u ∈ Rk .
such that u 2 = 1. and uT bj = 0. for all j = 1, · · · , k − 1..
10.10 Given A ∈ Rm×n ., let a1 , · · · , am ∈ Rn . be the column vectors of AT .. Prove
that
m
. AT 2F A 2
F = aj 2
2.
j =1
Appendix A
Related Topics

In this appendix, we review some related topics and results from calculus and linear
algebra.

A.1 Calculus

We say an infinite sequence {xn }∞n=1 ⊂ R. converges to c if for any ε > 0., there
exists N ∈ N. such that |xn − c| < ε. for all n > N .. In this case, we denote

. lim xn = c.
n→∞

The sequence {xn }∞n=1 ⊂ R. is called a Cauchy sequence if for any ε > 0., there
exists N ∈ N. such that |xn − xm | < ε. for all n > N . and m > N ..
Theorem A.1 The real line R. is complete; namely, the following statements are
equivalently valid:
(i) Every Cauchy sequence converges.
(ii) Every bounded and monotone sequence is convergent.
(iii) Every bounded sequence has a convergent subsequence.
(iv) Every non-empty bounded set has a supremum (least upper bound) and an
infimum (greatest lower bound)
Given a non-empty set A ⊂ R., we say a ∈ R. is a limit point of A if there exists an
infinite sequence {xn }∞n=1 ⊂ A. that converges to a: xn → a . as n → ∞.. The closure
of A, denoted by Ā., is the set containing all limit points of A. Since any point in A
is a limit point (by choosing xn = a . for all n = 1, 2, · · · .), we have A ⊂ Ā.. We say
A is closed if A = Ā..

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 103
X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
104 A Related Topics

Theorem A.2 A set A ⊂ R. is a closed and bounded set if and only if it is compact;
namely, any sequence in A has a subsequence that converges to a point in A.
Let f be a (real-valued) function defined on an open interval I = (a, b) ⊂ R.,
where a could be either finite or − ∞., and b could be either finite or ∞.. Given
c ∈ I ., we say f (x) → L. as x → c. if for any small ε > 0., there exists δ > 0. such
that 0 < |x − c| < δ . implies |f (x) − f (c)| < ε .. The value L is called the limit
of f (x). at c. The one-sided limits of f (x). at c can be defined in a similar manner.
We say f (x) → L. as x → c+ . (respectively, x → c− .) if for any small ε > 0., there
exists δ > 0. such that 0 < x − c < δ . (respectively, − δ < x − c < 0.) implies
|f (x) − f (c)| < ε .. We use

. lim f (x), lim f (x), lim f (x),

x→c x→c+ x→c−

respectively, to denote the limit, right limit, and left limit of f (x). at c. It is evident
that the limit exists if and only if both one-sided limits exist and equal. We say f is
continuous at c if

. lim f (x) = f (c).

x→c

One can prove that f is continuous at c if and only if for any convergent sequence
{xn }∞
n=1 ⊂ R. with limn→∞ xn = c., we have limn→∞ f (xn ) = f (c).. We say f is
differentiable at c if the limit

f (x) − f (c)
. lim
x→c x−c

exists, and in this case, we denote the value of the limit by f (c). and call it the
derivative of f at c. If f is continuous for each x ∈ I ., then we say f is continuous
on I . The set of all continuous functions on I is denoted by C(I ) = C(I, R).. When
the end points a and b are finite, we say f is continuous on I¯ = [a, b]. if f is
continuous on I = (a, b). and

. lim f (x) = f (a), lim f (x) = f (b).

x→a + x→b−

We say f is differentiable on I if f is differentiable at each x ∈ i. and f ∈

C(I ).; in this case we say f ∈ C 1 (I ) = C 1 (I, R).. The higher order derivatives are
defined as derivatives of derivatives. We say f ∈ C n (I ). if the derivative functions
f , f , · · · , f (n) . are continuous on I .
Theorem A.3 Let I = (a, b). with a < b., and assume f ∈ C(I¯).. We have the
following results:
(i) The set f (I¯) := {f (x) : x ∈ I¯}. is compact (i.e., bounded and closed).
A Related Topics 105

(ii) (Extreme Value Theorem) Let M and m be the supremum and infimum of f (I¯).,
respectively. There exist xM , xm ∈ I¯. such that f (xM ) = M . and f (xm ) = m..
(iii) (Intermediate Value Theorem) For any r ∈ (m, M)., where M and m are given
as above, there exists c ∈ (a, b). such that f (c) = r ..
(iv) (Mean Value Theorem) If further f ∈ C 1 (I )., then there exists ξ ∈ (a, b). such
that f (ξ ) = [f (b) − f (a)]/(b − a)..
Differentiation is a linear operator from C 1 (I ). to C(I ).; namely

(c1 f1 + c2 f2 ) = c1 f1 + c2 f2 ,
. (A.1)

for any c1 , c2 ∈ R. and f1 , f2 ∈ C 1 (I ).. The product rule and the quotient rule of
differentiation are

(f g) = f g + f g , (f/g) = (f g − fg )/g 2 ,
. (A.2)

for f, g ∈ C 1 (I ).. If g ∈ C 1 (I ). and f ∈ C 1 (g(I ))., we also have the chain rule

[f (g(x))] = f (g(x))g (x).

. (A.3)

Given f ∈ C(I )., we say F ∈ C 1 (I ). is an antiderivative of f if F (x) = f (x). for

all x ∈ (a, b).. Note that if F (x). is an antiderivative of f (x)., so is F (x) + C . for any
constant C. If further F ∈ C(I¯)., then we have the fundamental theorem of calculus

b b
. f (x)dx = F (x) = F (b) − F (a), (A.4)
a a

where the left-hand side is called a definite integral of f (x). over the interval I =
(a, b).. If f, g ∈ C 1 (I ) ∩ C(I¯)., then the product rule of differentiation (A.2) implies
the useful formula of integration by parts:

b b b
. f (x)g(x)dx = f (x)g(x) − f (x)g (x)dx. (A.5)
a a a

Theorem A.4 (Taylor’s Theorem) Let f ∈ C (n+1) (I ). with I = (a, b).. For any
x, x0 ∈ (a, b)., we have f (x) = Pn (x) + Rn (x)., where
n
f (k) (x0 )
Pn (x) =
. (x − x0 )k (A.6)
k!
k=0

is the Taylor polynomial, and

x (x − t)n (n+1)
Rn (x) =
. f (t)dt (A.7)
x0 n!
106 A Related Topics

is the remainder/error term. Moreover, there exists ξ ∈ (a, b). (depending on x and
x0 .) such that

(x − x0 )n+1 (n+1)
Rn (x) =
. f (ξ ). (A.8)
(n + 1)!

A.2 Linear Algebra

A matrix A ∈ Cm×n . is a set of m×n. complex numbers with the following structure:
⎛ ⎞
A11 · · · A1n
⎜ . .. ⎟ .
.A = ⎝ .
. . ⎠
Am1 · · · Amn

The Hermitian of a matrix A ∈ Cm×n . is denoted as AH ∈ Cn×m . such that

AH
ij = Aj i . (namely, AHij . is the complex conjugate of Aj i .) for all i = 1, · · · , n. and
j = 1, · · · , m.. A square matrix A ∈ Cn×n . is Hermitian if A = AH .. The transpose
of a real matrix A ∈ Rm×n . is denoted as AT ∈ Rn×m . such that ATij = Aj i . for all
i = 1, · · · , n. and j = 1, · · · , m.. A real square matrix A ∈ Rn×n . is symmetric
if A = AT .. The matrix multiplication of A ∈ Cm×n . and B ∈ Cn×d . is defined as
AB ∈ Cm×d . such that
n
(AB)ij =
. Aik Bkj ,
k=1

for all i = 1, · · · , m. and j = 1, · · · , d .. An identity matrix I ∈ Cn×n . is a special

matrix defined as

1, i = j,
Iij = δij =
.
0, i j,

for i, j = 1, · · · , n.. It is easily seen that AI = A. and I B = B . for any A ∈ Cm×n .

and B ∈ Cn×d .. We say a square matrix A ∈ Cn×n . is invertible/nonsingular if
there exists B ∈ Cn×n . such that AB = BA = I .. A square matrix A ∈ Cn×n . is
unitary if AAH = AH A = I .. Moreover, we say a square matrix A ∈ Cn×n . is
normal if AAH = AH A.. The diagonal, Hermitian, and unitary matrices are special
examples of normal matrices. The determinant of a square matrix A ∈ Cn×n . is
defined recursively as
A Related Topics 107

A11 , n = 1,
. det(A) = n j −1 A det(A(1,j ) ), n > 1,
j =1 (−1) 1j

where A(i,j ) ∈ C(n−1)×(n−1) . is the submatrix obtained by deleting the i-th row and
the j -th column in A. The trace of a square matrix A ∈ Cn×n . is defined as the sum
of diagonal terms:
n
. tr(A) := Ajj .
j =1

For any A, B ∈ Cn×n ., it follows from the definition that

n n
. tr(AB) = tr(BA) = Aj k Bkj .
j =1 k=1

A (column) vector v ∈ Cn = Cn×1 . is a special matrix. We say a set of

vectors v1 , · · · , vn ∈ Cn . is linearly dependent if and only if there exist constants
c1 , · · · , cn ∈ C., not identically zero, such that

.c1 v1 + · · · + cn vn = 0 ∈ Cn .

In other words, the vectors v1 , · · · , vn . are linearly independent if and only if the
above equation has only one trivial solution c1 = · · · = cn = 0.. If the vectors
v1 , · · · , vk ∈ Cn . are linearly independent, then we can define a k-dimensional
subspace Vk ∈ Cn . spanned by these vectors as follows:

Vk = span{v1 , · · · , vk } := {c1 v1 + · · · + ck vk , c1 , · · · , ck ∈ C}.

.

In particular, if k = n. and v1 , · · · , vn ∈ Cn . are linearly independent, then Vn =

span{v1 , · · · , vn } = Cn .. We call v1 , · · · , vk . a basis of Vk . in the sense that any
vector v ∈ Vk . can be expressed as a linear combination of the basis. If these vectors
are orthogonal (i.e., vjH vl = 0. for 1 ≤ j < l ≤ k .), then the basis is orthogonal.
Moreover, we say the basis is orthonormal if in addition to orthogonality, the vectors
are normalized as vj 2 = 1. for j = 1, · · · , k ..
Theorem A.5 Let A ∈ Cn×n .. The following statements are equivalent:
(i) A is invertible (nonsingular).
(ii) det(A) 0..
(iii) The linear homogeneous equation Ax = 0 ∈ Cn . does not have any nontrivial
solution x 0. in Cn ..
(iv) The linear inhomogeneous equation Ax = b. has a unique solution x ∈ Cn . for
any given b ∈ Cn ..
108 A Related Topics

(v) The columns of A, considering as a set of n column vectors, are linearly

independent in C n ..
The characteristic polynomial of a square matrix A ∈ Cn×n . is pn (λ) = det(λI −
A) = λn + · · · ., whose zeros are called the eigenvalues of A. According to the
fundamental theorem of algebra, there are n eigenvalues (counting multiplicity),
denoted by λ1 , · · · , λn .. The spectral radius of A is the largest modulus of the
eigenvalues:

ρ(A) = max |λj |.

.
j

For each eigenvalue λ. of A, there exists a nonzero eigenvector x ∈ Cn . such

that Ax = λx .. If A ∈ Cn×n . is Hermitian, then (x H Ax)H = x H Ax ∈ R. for any
x ∈ Cn .. Hence, the eigenvalues of a Hermitian matrix are real. We say a Hermitian
matrix A ∈ Cn×n . is positive (nonnegative) definite if x H Ax > 0. (x H Ax ≥ 0.) for all
x ∈ Cn . with x 0.. A Hermitian matrix A ∈ Cn×n . is positive (nonnegative) definite
if and only if all eigenvalues of A are positive (nonnegative). A nonnegative definite
Hermitian matrix is also called positive semidefinite in the literature. A Hermitian
matrix A ∈ Cn×n . is nonnegative definite if and only if there exists B ∈ Cn×n . such
that A = B H B .. If further A is positive definite, then B is nonsingular/invertible.
. on C . is a function satisfying the following three conditions:
A vector norm n

(i) (Positivity) u > 0. for all u ∈ Cn . with u 0..

(ii) (Homogeneity) cu c u . for all c ∈ C . and u ∈ Cn ..
(iii) (Triangle inequality) u + v u v . for all u, v ∈ Cn ..
The most commonly used vector norms are l p .-norm with p = 1, 2, ∞.:

n n
. u 1 := |ui |, u 2 := |ui |2 , u ∞ := max |ui |.
i
i=1 i=1

An inner product (·, ·). is a bilinear function on Cn × Cn . satisfying the following

three conditions:
(i) (Positivity) (u, u) > 0. for all u ∈ Cn . with u 0..
(ii) (Symmetry) (u, v) = (v, u). for all u, v ∈ Cn ..
(iii) (Linearity) (au+bv, w) = a(u, w)+b(v, w). for all a, b ∈ C. and u, v, w ∈ Cn ..
√
It can be verified that u (u, u). is a norm on Cn ..
β on C . are equivalent if there exist positive constants
Two norms . and .
n
α
C1 > 0. and C2 > 0. such that C1 u α u β ≤ C2 u α . for all u ∈ Cn ..

√ on C . are equivalent and continuous. In particular, we

Theorem A.6 All norms n

have u ∞
√ u 2 ≤ n u ∞ ., u ∞ u 1 ≤ n u ∞ ., and u 2 u 1 ≤
n u 2 . for all u ∈ Cn ..
A Related Topics 109

A matrix norm . on Cm×n . is a function satisfying the following three

conditions:
(i) (Positivity) A > 0. for all A ∈ Cm×n . with A 0..
(ii) (Homogeneity) cA c A . for all c ∈ C . and A ∈ Cm×n ..
(iii) (Triangle inequality) A + B A B . for all A, B ∈ Cm×n ..
A typical matrix norm is the Frobenius norm, which is defined by treating a matrix
in Cm×n . as a vector of dimension mn and then applying the l 2 .-norm in the mn-
dimensional vector space:

m n
. A F := |Aj k |2 = tr(AH A) = tr(AAH ), A ∈ Cm×n .
j =1 k=1

In many cases, we mainly focus on matrix norms of square matrices. In this case,
we always scale the norm such that it satisfies

. AB A B , (A.9)

for all A, B ∈ Cn×n .. Let λ. be an eigenvalue of A ∈ Cn×n . with a nonzero

eigenvector x ∈ Cn .. Denote B = (x, · · · , x) ∈ Cn×n .. We have

|λ B
. λB AB A B ,

which implies |λ A .. Since λ. is an arbitrary eigenvalue, we have ρ(A) A ..

This implies that the spectral radius is a lower bound of the (scaled) matrix norm. On
the other hand, for any A ∈ Cn×n . and ε > 0., there exists a matrix norm (depending
on A and ε .) such that A < ρ(A) + ε..
. on C ., there exists a natural/induced matrix norm
Given any vector norm n

. on C
n×n . that is defined by

Au
. A sup = sup Av ,
u∈Cn , u 0 u v∈Cn , v 1

for any A ∈ C n×n .. The most commonly used matrix norms are l p .-norms with
p = 1, 2, ∞.:
n
. A 1 := sup Av 1 = max |Aj k |,
v∈Cn , v 1 =1 1≤k≤n
j =1

A 2 := sup Av 2 = ρ(AH A),

v∈Cn , v 2 =2
110 A Related Topics

n
A ∞ := sup Av ∞ = max |Aj k |.
v∈Cn , v ∞ =1
1≤j ≤n
k=1

Theorem A.7 Let A ∈ C n×n .. The following statements are equivalent:

(i) ρ(A) < 1..
(ii) limk→∞ Ak = 0..
(iii) limk→∞ Ak v = 0. for any v ∈ Cn ..
(iv) There exists a matrix norm . such that A < 1..

(v) The series ∞ k=0 A k . converges.

Appendix B
Hints to Selected Exercise Problems

We provide solution hints to some of the exercise problems.

1.6 The general solution to the second-order difference equation

Dn+1 = 2 cos θ Dn − Dn−1

.

is

Dn = aeiθ + be−iθ ,
.

where a and b are constants determined by the initial values.

1.7 If xj = xk . for some j k ., then det(B) = 0.. Hence, j k (xj −xk ). is a factor
of det(B).. On the other hand, det(B). is a polynomial in (x1 , · · · , xn ) ∈ Rn .
of degree n(n − 1)/2.. Therefore,

. det(B) = c (xj − xk ),
1≤k<j ≤n

where c is a constant.
2.4 Define f (t) := J (p+t (q−p)). for t ∈ [0, 1]., and then apply Taylor’s theorem.
3.7 Fix u ∈ I .. For any sufficiently small ε > 0., the function g(h) := [f (u + h) −
f (u)]/ h. is monotone and bounded for h ∈ (0, ε]. and also for h ∈ [−ε, 0)..
4.7 Apply Taylor’s theorem to the function f (α) := J (p − αJ (p))..
5.3 Note that (x1 ∨ x2 ) ∧ x3 = 1. if and only if x1 ∨ x2 = 1. and x3 = 1., which is
equivalent to x1 ∧ x3 = 1. or x2 ∧ x3 = 1.. The other identity can be proved in
a similar manner.
6.1 If i k ., then Cij . is independent of Akl .. If i = k ., then from

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 111
X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
112 B Hints to Selected Exercise Problems

m
Cij =
. Ail Blj
l=1

we obtain
∂Cij
. = Blj .
Ail

The other identity follows similarly.

√ √
9.11 If γ ≥ 2 β ., then |γ + γ 2 − 4β| = γ + γ 2 − 4β ≥ γ ≥ 2 β .. If
√ √
γ < 2 β ., then |γ + γ 2 − 4β| = 2 β ..
10.3 Consider the Schur decomposition A = U T U H ., where U is unitary and T is
upper triangular. Since A is normal, we have T T H = T H T .. One can verify
that all off-diagonal terms of T vanish.
Appendix C
Further Readings

Machine learning has strong connections to statistics, optimization, numerical

analysis, and mathematical analysis. In this book, we briefly discuss some selected
topics and provide pointers for further study.
Interested readers are referred to more detailed studies on regression [5, 7],
variance inflation factor [13], convex optimization [1, 4], LASSO [16], neural
networks [3], deep learning [9], batch normalization [12], support vector machines
[2, 18], the conjugate gradient method [10], the heavy ball method [17], singular
value decomposition [6, 8], and principal component analysis [11, 14, 15].
The field of machine learning is rapidly evolving, and the references cited in
the above paragraph are by no means comprehensive. Due to the limits of our
knowledge, we only provide a selection of relevant references and apologize for
omitting many other important works.

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 113
X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
References

1. S. Boyd and L. Vandenberghe, Convex Optimization. Cambridge University Press, Cambridge,

UK, 2004.
2. C. Cortes and V. Vapnik, Support-vector networks, Mach. Learn. 20 (1995), 273–297.
3. G. Cybenko, Approximation by superpositions of a sigmoidal function, Math. Control Signals
Syst. 2 (1989), 303–314.
4. I. Daubechies, M. Defrise, and C. De Mol, An iterative thresholding algorithm for linear inverse
problems with a sparsity constraint, Commun. Pure Appl. Math. 57 (2004), 1413–1457.
5. N. R. Draper and H. Smith, Applied Regression Analysis, 3rd ed. Wiley, New York, 1998.
6. C. Eckart and G. Young, The approximation of one matrix by another of lower rank,
Psychometrika 1 (1936), 211–218.
7. D. A. Freedman, Statistical Models: Theory and Practice. Cambridge University Press,
Cambridge, UK, 2009.
8. G. H. Golub and C. F. Van Loan, Matrix Computations, 4th ed. Johns Hopkins University Press,
Baltimore, MD, 2013.
9. I. Goodfellow, Y. Bengio, and A. Courville, Deep Learning. MIT Press, Cambridge, MA, 2016.
10. M. R. Hestenes and E. Stiefel, Methods of conjugate gradients for solving linear systems, J.
Res. Nat. Bur. Stand. 49 (1952), 409–436.
11. H. Hotelling, Analysis of a complex of statistical variables into principal components, J. Educ.
Psychol. 24 (1933), 417–441, 498–520.
12. S. Ioffe and C. Szegedy, Batch normalization: Accelerating deep network training by reducing
internal covariate shift, in Proc. 32nd Int. Conf. Mach. Learn. (ICML), 2015, pp. 448–456.
13. G. James, D. Witten, T. Hastie, and R. Tibshirani, An Introduction to Statistical Learning, 2nd
ed. Springer, New York, 2021.
14. I. T. Jolliffe, Principal Component Analysis, 2nd ed. Springer, New York, 2002.
15. I. T. Jolliffe and J. Cadima, Principal component analysis: A review and recent developments,
Phil. Trans. R. Soc. A 374 (2016), 20150202.
16. R. Tibshirani, Regression shrinkage and selection via the lasso, J. R. Stat. Soc. Ser. B 58 (1996),
267–288.
17. B. T. Polyak, Some methods of speeding up the convergence of iteration methods, USSR
Comput. Math. Math. Phys. 4 (1964), 1–17.
18. V. N. Vapnik, Statistical Learning Theory. Wiley, New York, 1998.

© The Author(s), under exclusive license to Springer Nature Switzerland AG 2026 115
X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
Index

B cost, 48
Backward propagation, 64, 71 differentiable, 104
Batch normalization, 69 Fréchet differentiable, 12
gradient, 12
inner product, 83, 108
C likelihood, 19
Cauchy sequence, 103 multilinear, 3
Coefficient of determination, 24, 32 normalized, 3
Conjugate gradient method, 83 sigmoid, 46
Convex optimization, 34
Covariance, 19
G
Gradient descent method, 49, 58
D Gram-Schmidt process, 84, 93
Data
dimension, 17
nontrivial, 76 H
separable, 76 Hyperplane, 76
size, 17

I
E Input, 17
Error, 19 binary, 45
Iterative method, 39
heavy ball method, 88
F Jacobi method, 87
Function one-step, 86
activation, 60, 69 two-step, 88
antisymmetric, 3
composite, 57
continuous, 104 K
convex, 34 Kronecker delta, 77

© The Editor(s) (if applicable) and The Author(s), under exclusive license to 117
Springer Nature Switzerland AG 2026
X.-S. Wang, C. Wang, Machine Learning in Data Processing, Forum for
Interdisciplinary Mathematics, [Link]
118 Index

L l 2 ., 99
Layer l p ., 108
hidden, 60, 63 Frobenius, 99, 109
input, 60, 63 matrix, 109
output, 60, 63 vector, 75, 108
Least absolute shrinkage and selection operator Normal distribution, 19, 31
(LASSO), 34 multivariate, 19
Least squares approximation, 20
Limit point, 103
Linear correlation, 25 O
Linear system, 27, 83, 86 Objective function, 20
Logic gate, 45 Observation, 18, 75
Optimization problem, 20, 48, 76
constrained, 77
M minimax, 87, 88
Matrix, 1 Output, 17
characteristic polynomial, 108 binary, 45
conjugate transpose, 93 Overfitting, 32
determinant, 3, 106
eigenvalue, 108
eigenvector, 108 P
factorization, 94 Parameter, 18, 75
Hermitian, 106 feasible, 76
Hessian, 13 learning rate, 39, 50
identity, 2, 106 normalized feasible, 76
invertible/nonsingular, 2, 106 regularization, 41
Jacobian, 13 tuning/hyperparameter, 33
minor, 1 Prediction, 18, 46
multiplication, 106 Principal component analysis (PCA), 95
nonnegative definite/positive semidefinite, Projection, 21, 84, 96
108 Pythagorean theorem, 21, 97
normal, 106
orthogonal, 95
positive definite, 83, 108
principal component, 96 R
product, 2 Rayleigh quotient, 95
rank, 9 Regression
singular value, 95 linear, 17
spectral radius, 108 nonlinear, 48
symmetric, 2, 106 Ridge, 32
trace, 2, 107 Regularization, 31
transpose, 2, 106 l 1 ., 40
unitary, 93, 106 l 2 ., 32, 40
Maximum likelihood estimation, 19 Residual, 84
Mean, 19, 31

S
N Sample variance, 24
Neural network Schur decomposition, 93
deep, 63 Set
depth, 63 closed, 103
shallow, 60 compact, 104
Norm feasible, 75
l 1 ., 34 normalized feasible, 76
Index 119

Singular value decomposition (SVD), 95 V

Steepest descent method, 39 Variance, 19, 31
Subderivative, 36 Variance inflation factor, 25, 32
Subdifferential, 36 Vectors, 1
Sum of squared errors (SSE), 19 basis, 107
Support vector machine (SVM), 75 linearly dependent, 11, 107
linearly independent, 11, 107
orthogonal basis, 107
T orthonormal, 96
Taylor expansion, 13, 21, 33 orthonormal basis, 99, 107
Taylor’s theorem, 14, 105 standard basis, 2