1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
|
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.spark.ml.optim
import org.apache.spark.internal.Logging
import org.apache.spark.ml.feature.Instance
import org.apache.spark.ml.linalg._
import org.apache.spark.rdd.RDD
/**
* Model fitted by [[IterativelyReweightedLeastSquares]].
* @param coefficients model coefficients
* @param intercept model intercept
* @param diagInvAtWA diagonal of matrix (A^T * W * A)^-1 in the last iteration
* @param numIterations number of iterations
*/
private[ml] class IterativelyReweightedLeastSquaresModel(
val coefficients: DenseVector,
val intercept: Double,
val diagInvAtWA: DenseVector,
val numIterations: Int) extends Serializable
/**
* Implements the method of iteratively reweighted least squares (IRLS) which is used to solve
* certain optimization problems by an iterative method. In each step of the iterations, it
* involves solving a weighted least squares (WLS) problem by [[WeightedLeastSquares]].
* It can be used to find maximum likelihood estimates of a generalized linear model (GLM),
* find M-estimator in robust regression and other optimization problems.
*
* @param initialModel the initial guess model.
* @param reweightFunc the reweight function which is used to update offsets and weights
* at each iteration.
* @param fitIntercept whether to fit intercept.
* @param regParam L2 regularization parameter used by WLS.
* @param maxIter maximum number of iterations.
* @param tol the convergence tolerance.
*
* @see [[http://www.jstor.org/stable/2345503 P. J. Green, Iteratively Reweighted Least Squares
* for Maximum Likelihood Estimation, and some Robust and Resistant Alternatives,
* Journal of the Royal Statistical Society. Series B, 1984.]]
*/
private[ml] class IterativelyReweightedLeastSquares(
val initialModel: WeightedLeastSquaresModel,
val reweightFunc: (Instance, WeightedLeastSquaresModel) => (Double, Double),
val fitIntercept: Boolean,
val regParam: Double,
val maxIter: Int,
val tol: Double) extends Logging with Serializable {
def fit(instances: RDD[Instance]): IterativelyReweightedLeastSquaresModel = {
var converged = false
var iter = 0
var model: WeightedLeastSquaresModel = initialModel
var oldModel: WeightedLeastSquaresModel = null
while (iter < maxIter && !converged) {
oldModel = model
// Update offsets and weights using reweightFunc
val newInstances = instances.map { instance =>
val (newOffset, newWeight) = reweightFunc(instance, oldModel)
Instance(newOffset, newWeight, instance.features)
}
// Estimate new model
model = new WeightedLeastSquares(fitIntercept, regParam, standardizeFeatures = false,
standardizeLabel = false).fit(newInstances)
// Check convergence
val oldCoefficients = oldModel.coefficients
val coefficients = model.coefficients
BLAS.axpy(-1.0, coefficients, oldCoefficients)
val maxTolOfCoefficients = oldCoefficients.toArray.reduce { (x, y) =>
math.max(math.abs(x), math.abs(y))
}
val maxTol = math.max(maxTolOfCoefficients, math.abs(oldModel.intercept - model.intercept))
if (maxTol < tol) {
converged = true
logInfo(s"IRLS converged in $iter iterations.")
}
logInfo(s"Iteration $iter : relative tolerance = $maxTol")
iter = iter + 1
if (iter == maxIter) {
logInfo(s"IRLS reached the max number of iterations: $maxIter.")
}
}
new IterativelyReweightedLeastSquaresModel(
model.coefficients, model.intercept, model.diagInvAtWA, iter)
}
}
|
