Coverage for picos/expressions/exp

1# ------------------------------------------------------------------------------

3# Based on the original picos.expressions module by Guillaume Sagnol.

5# This file is part of PICOS.

7# PICOS is free software: you can redistribute it and/or modify it under the

8# terms of the GNU General Public License as published by the Free Software

9# Foundation, either version 3 of the License, or (at your option) any later

10# version.

11#

12# PICOS is distributed in the hope that it will be useful, but WITHOUT ANY

13# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR

14# A PARTICULAR PURPOSE. See the GNU General Public License for more details.

15#

16# You should have received a copy of the GNU General Public License along with

17# this program. If not, see <http://www.gnu.org/licenses/>.

18# ------------------------------------------------------------------------------

20"""Implements :class:`Entropy` and :class:`NegativeEntropy`."""

22# TODO: Common base class for Entropy and NegativeEntropy.

24import math

25import operator

26from collections import namedtuple

28import cvxopt

29import numpy

31from .. import glyphs

32from ..apidoc import api_end, api_start

33from ..caching import cached_selfinverse_unary_operator, cached_unary_operator

34from ..constraints import KullbackLeiblerConstraint

35from .data import convert_and_refine_arguments, convert_operands, cvx2np

36from .exp_affine import AffineExpression

37from .expression import Expression, refine_operands, validate_prediction

39_API_START = api_start(globals())

40# -------------------------------

43class Entropy(Expression):

44 r"""Entropy or negative relative entropy of an affine expression.

46 Negative relative entropy is also known as the perspective of the logarithm.

48 :Definition:

50 Let :math:`x` be an :math:`n`-dimensional real affine expression.

52 1. If no additional expression :math:`y` is given, this is the entropy

54 .. math::

56 -\sum_{i = 1}^n \operatorname{vec}(x)_i

57 \log(\operatorname{vec}(x)_i).

59 2. If an additional affine expression :math:`y` of same shape as :math:`x`

60 is given, this is the negative relative entropy (or logarithmic

61 perspective)

63 .. math::

65 -\sum_{i = 1}^n \operatorname{vec}(x)_i

66 \log\left(

67 \frac{\operatorname{vec}(x)_i}{\operatorname{vec}(y)_i}

68 \right)

69 &= -\sum_{i = 1}^n \operatorname{vec}(x)_i

70 \left[

71 \log(\operatorname{vec}(x)_i) - \log(\operatorname{vec}(y)_i)

72 \right] \\

73 &= \sum_{i = 1}^n \operatorname{vec}(x)_i

74 \left[

75 \log(\operatorname{vec}(y)_i) - \log(\operatorname{vec}(x)_i

76 \right] \\

77 &= \sum_{i = 1}^n \operatorname{vec}(x)_i

78 \log\left(

79 \frac{\operatorname{vec}(y)_i}{\operatorname{vec}(x)_i}

80 \right).

82 .. warning::

84 When you pose a lower bound on this expression, then PICOS enforces

85 :math:`x \geq 0` through an auxiliary constraint during solution search.

86 When an additional expression :math:`y` is given, PICOS enforces

87 :math:`y \geq 0` as well.

88 """

90 # --------------------------------------------------------------------------

91 # Initialization and factory methods.

92 # --------------------------------------------------------------------------

94 @convert_and_refine_arguments("x", "y", allowNone=True)

95 def __init__(self, x, y=None):

96 """Construct an :class:`Entropy`.

98 :param x: The affine expression :math:`x`.

99 :type x: ~picos.expressions.AffineExpression

100 :param y: An additional affine expression :math:`y`. If necessary, PICOS

101 will attempt to reshape or broadcast it to the shape of :math:`x`.

102 :type y: ~picos.expressions.AffineExpression

103 """

104 if not isinstance(x, AffineExpression):

105 raise TypeError("Can only take the elementwise logarithm of a real "

106 "affine expression, not of {}.".format(type(x).__name__))

107

108 if y is not None:

109 if not isinstance(y, AffineExpression):

110 raise TypeError("The additional parameter y must be a real "

111 "affine expression, not {}.".format(type(x).__name__))

112 elif x.shape != y.shape:

113 y = y.reshaped_or_broadcasted(x.shape)

114

115 if y.is1:

116 y = None

117

118 self._x = x

119 self._y = y

120

121 if y is None:

122 typeStr = "Entropy"

123 if len(x) == 1:

124 symbStr = glyphs.neg(glyphs.mul(x.string, glyphs.log(x.string)))

125 else:

126 symbStr = glyphs.neg(

127 glyphs.dotp(x.string, glyphs.log(x.string)))

128 else:

129 typeStr = "Logarithmic Perspective"

130 if len(x) == 1:

131 symbStr = glyphs.mul(

132 x.string, glyphs.log(glyphs.div(y.string, x.string)))

133 else:

134 symbStr = glyphs.dotp(x.string,

135 glyphs.sub(glyphs.log(y.string), glyphs.log(x.string)))

136

137 Expression.__init__(self, typeStr, symbStr)

138

139 # --------------------------------------------------------------------------

140 # Abstract method implementations and method overridings, except _predict.

141 # --------------------------------------------------------------------------

142

143 def _get_refined(self):

144 if self._x.constant and (self._y is None or self._y.constant):

145 return AffineExpression.from_constant(self.value, 1, self._symbStr)

146 else:

147 return self

148

149 Subtype = namedtuple("Subtype", ("argdim", "y"))

150

151 def _get_subtype(self):

152 return self.Subtype(len(self._x), self._y is not None)

153

154 def _get_value(self):

155 x = numpy.ravel(cvx2np(self._x._get_value()))

156

157 if self._y is None:

158 s = -numpy.dot(x, numpy.log(x))

159 else:

160 y = numpy.ravel(cvx2np(self._y._get_value()))

161 s = numpy.dot(x, numpy.log(y / x))

162

163 return cvxopt.matrix(s)

164

165 @cached_unary_operator

166 def _get_mutables(self):

167 if self._y is None:

168 return self._x._get_mutables()

169 else:

170 return self._x._get_mutables().union(self._y.mutables)

171

172 def _is_convex(self):

173 return False

174

175 def _is_concave(self):

176 return True

177

178 def _replace_mutables(self, mapping):

179 return self.__class__(self._x._replace_mutables(mapping),

180 None if self._y is None else self._y._replace_mutables(mapping))

181

182 def _freeze_mutables(self, freeze):

183 return self.__class__(self._x._freeze_mutables(freeze),

184 None if self._y is None else self._y._freeze_mutables(freeze))

185

186 # --------------------------------------------------------------------------

187 # Python special method implementations, except constraint-creating ones.

188 # --------------------------------------------------------------------------

189

190 @cached_selfinverse_unary_operator

191 def __neg__(self):

192 return NegativeEntropy(self._x, self._y)

193

194 @convert_operands(scalarRHS=True)

195 @refine_operands()

196 def __add__(self, other):

197 other_str = other.string

198

199 if isinstance(other, AffineExpression):

200 if other.is0:

201 return self

202

203 other = Entropy(other, other * math.e)

204

205 if isinstance(other, Entropy):

206 if self._y is None and other._y is None:

207 entropy = Entropy(self._x.vec // other._x.vec)

208 elif self._y is not None and other._y is None:

209 one = AffineExpression.from_constant(1.0, (other.n, 1))

210 entropy = Entropy(

211 self._x.vec // other._x.vec, self._y.vec // one)

212 elif self._y is None and other._y is not None:

213 one = AffineExpression.from_constant(1.0, (self.n, 1))

214 entropy = Entropy(

215 self._x.vec // other._x.vec, one // other._y.vec)

216 else:

217 entropy = Entropy(

218 self._x.vec // other._x.vec, self._y.vec // other._y.vec)

219

220 entropy._symbStr = glyphs.clever_add(self.string, other_str)

221

222 return entropy

223

224 return Expression.__add__(self, other)

225

226 @convert_operands(scalarRHS=True)

227 @refine_operands()

228 def __radd__(self, other):

229 if isinstance(other, (AffineExpression, Entropy)):

230 entropy = Entropy.__add__(self, other)

231

232 if entropy is not self:

233 entropy._symbStr = glyphs.clever_add(other.string, self.string)

234

235 return entropy

236

237 return Expression.__radd__(self, other)

238

239 @convert_operands(scalarRHS=True)

240 @refine_operands()

241 def __sub__(self, other):

242 if isinstance(other, (AffineExpression, NegativeEntropy)):

243 return self + (-other)

244

245 return Expression.__sub__(self, other)

246

247 # --------------------------------------------------------------------------

248 # Methods and properties that return expressions.

249 # --------------------------------------------------------------------------

250

251 @property

252 def x(self):

253 """The expression :math:`x`."""

254 return self._x

255

256 @property

257 def y(self):

258 """The additional expression :math:`y`, or :obj:`None`."""

259 return self._y

260

261 # --------------------------------------------------------------------------

262 # Methods and properties that describe the expression.

263 # --------------------------------------------------------------------------

264

265 @property

266 def n(self):

267 """Length of :attr:`x`."""

268 return len(self._x)

269

270 # --------------------------------------------------------------------------

271 # Constraint-creating operators, and _predict.

272 # --------------------------------------------------------------------------

273

274 @classmethod

275 def _predict(cls, subtype, relation, other):

276 assert isinstance(subtype, cls.Subtype)

277

278 if relation == operator.__ge__:

279 if issubclass(other.clstype, AffineExpression) \

280 and other.subtype.dim == 1:

281 return KullbackLeiblerConstraint.make_type(

282 argdim=subtype.argdim)

283

284 return NotImplemented

285

286 @convert_operands(scalarRHS=True)

287 @validate_prediction

288 @refine_operands()

289 def __ge__(self, other):

290 if isinstance(other, AffineExpression):

291 return KullbackLeiblerConstraint(-self, -other)

292 else:

293 return NotImplemented

294

295

296class NegativeEntropy(Expression):

297 r"""Negative or relative entropy of an affine expression.

298

299 Relative entropy is also known as the Kullback-Leibler divergence.

300

301 :Definition:

302

303 Let :math:`x` be an :math:`n`-dimensional real affine expression.

304

305 1. If no additional expression :math:`y` is given, this is the negative

306 entropy

307

308 .. math::

309

310 \sum_{i = 1}^n \operatorname{vec}(x)_i

311 \log(\operatorname{vec}(x)_i).

312

313 2. If an additional affine expression :math:`y` of same shape as :math:`x`

314 is given, this is the relative entropy (or Kullback-Leibler divergence)

315

316 .. math::

317

318 \sum_{i = 1}^n \operatorname{vec}(x)_i

319 \log\left(

320 \frac{\operatorname{vec}(x)_i}{\operatorname{vec}(y)_i}

321 \right)

322 = \sum_{i = 1}^n \operatorname{vec}(x)_i

323 \left[

324 \log(\operatorname{vec}(x)_i) - \log(\operatorname{vec}(y)_i)

325 \right].

326

327 .. warning::

328

329 When you pose an upper bound on this expression, then PICOS enforces

330 :math:`x \geq 0` through an auxiliary constraint during solution search.

331 When an additional expression :math:`y` is given, PICOS enforces

332 :math:`y \geq 0` as well.

333 """

334

335 # --------------------------------------------------------------------------

336 # Initialization and factory methods.

337 # --------------------------------------------------------------------------

338

339 @convert_and_refine_arguments("x", "y", allowNone=True)

340 def __init__(self, x, y=None):

341 """Construct a :class:`NegativeEntropy`.

342

343 :param x: The affine expression :math:`x`.

344 :type x: ~picos.expressions.AffineExpression

345 :param y: An additional affine expression :math:`y`. If necessary, PICOS

346 will attempt to reshape or broadcast it to the shape of :math:`x`.

347 :type y: ~picos.expressions.AffineExpression

348 """

349 if not isinstance(x, AffineExpression):

350 raise TypeError("Can only take the elementwise logarithm of a real "

351 "affine expression, not of {}.".format(type(x).__name__))

352

353 if y is not None:

354 if not isinstance(y, AffineExpression):

355 raise TypeError("The additional parameter y must be a real "

356 "affine expression, not {}.".format(type(x).__name__))

357 elif x.shape != y.shape:

358 y = y.reshaped_or_broadcasted(x.shape)

359

360 if y.is1:

361 y = None

362

363 self._x = x

364 self._y = y

365

366 if y is None:

367 typeStr = "Negative Entropy"

368 if len(x) == 1:

369 symbStr = glyphs.mul(x.string, glyphs.log(x.string))

370 else:

371 symbStr = glyphs.dotp(x.string, glyphs.log(x.string))

372 else:

373 typeStr = "Relative Entropy"

374 if len(x) == 1:

375 symbStr = glyphs.mul(

376 x.string, glyphs.log(glyphs.div(x.string, y.string)))

377 else:

378 symbStr = glyphs.dotp(x.string,

379 glyphs.sub(glyphs.log(x.string), glyphs.log(y.string)))

380

381 Expression.__init__(self, typeStr, symbStr)

382

383 # --------------------------------------------------------------------------

384 # Abstract method implementations and method overridings, except _predict.

385 # --------------------------------------------------------------------------

386

387 def _get_refined(self):

388 if self._x.constant and (self._y is None or self._y.constant):

389 return AffineExpression.from_constant(self.value, 1, self._symbStr)

390 else:

391 return self

392

393 Subtype = namedtuple("Subtype", ("argdim", "y"))

394

395 def _get_subtype(self):

396 return self.Subtype(len(self._x), self._y is not None)

397

398 def _get_value(self):

399 x = numpy.ravel(cvx2np(self._x._get_value()))

400

401 if self._y is None:

402 s = numpy.dot(x, numpy.log(x))

403 else:

404 y = numpy.ravel(cvx2np(self._y._get_value()))

405 s = numpy.dot(x, numpy.log(x / y))

406

407 return cvxopt.matrix(s)

408

409 @cached_unary_operator

410 def _get_mutables(self):

411 if self._y is None:

412 return self._x._get_mutables()

413 else:

414 return self._x._get_mutables().union(self._y.mutables)

415

416 def _is_convex(self):

417 return True

418

419 def _is_concave(self):

420 return False

421

422 def _replace_mutables(self, mapping):

423 return self.__class__(self._x._replace_mutables(mapping),

424 None if self._y is None else self._y._replace_mutables(mapping))

425

426 def _freeze_mutables(self, freeze):

427 return self.__class__(self._x._freeze_mutables(freeze),

428 None if self._y is None else self._y._freeze_mutables(freeze))

429

430 # --------------------------------------------------------------------------

431 # Python special method implementations, except constraint-creating ones.

432 # --------------------------------------------------------------------------

433

434 @cached_selfinverse_unary_operator

435 def __neg__(self):

436 return Entropy(self._x, self._y)

437

438 @convert_operands(scalarRHS=True)

439 @refine_operands()

440 def __add__(self, other):

441 other_str = other.string

442

443 if isinstance(other, AffineExpression):

444 if other.is0:

445 return self

446

447 other = NegativeEntropy(other, other / math.e)

448

449 if isinstance(other, NegativeEntropy):

450 if self._y is None and other._y is None:

451 negent = NegativeEntropy(self._x.vec // other._x.vec)

452 elif self._y is not None and other._y is None:

453 one = AffineExpression.from_constant(1.0, (other.n, 1))

454 negent = NegativeEntropy(

455 self._x.vec // other._x.vec, self._y.vec // one)

456 elif self._y is None and other._y is not None:

457 one = AffineExpression.from_constant(1.0, (self.n, 1))

458 negent = NegativeEntropy(

459 self._x.vec // other._x.vec, one // other._y.vec)

460 else:

461 negent = NegativeEntropy(

462 self._x.vec // other._x.vec, self._y.vec // other._y.vec)

463

464 negent._symbStr = glyphs.clever_add(self.string, other_str)

465

466 return negent

467

468 return Expression.__add__(self, other)

469

470 @convert_operands(scalarRHS=True)

471 @refine_operands()

472 def __radd__(self, other):

473 if isinstance(other, (AffineExpression, NegativeEntropy)):

474 negent = NegativeEntropy.__add__(self, other)

475

476 if negent is not self:

477 negent._symbStr = glyphs.clever_add(other.string, self.string)

478

479 return negent

480

481 return Expression.__radd__(self, other)

482

483 @convert_operands(scalarRHS=True)

484 @refine_operands()

485 def __sub__(self, other):

486 if isinstance(other, (AffineExpression, Entropy)):

487 return self + (-other)

488

489 return Expression.__sub__(self, other)

490

491 # --------------------------------------------------------------------------

492 # Methods and properties that return expressions.

493 # --------------------------------------------------------------------------

494

495 @property

496 def x(self):

497 """The expression :math:`x`."""

498 return self._x

499

500 @property

501 def y(self):

502 """The additional expression :math:`y`, or :obj:`None`."""

503 return self._y

504

505 # --------------------------------------------------------------------------

506 # Methods and properties that describe the expression.

507 # --------------------------------------------------------------------------

508

509 @property

510 def n(self):

511 """Length of :attr:`x`."""

512 return len(self._x)

513

514 # --------------------------------------------------------------------------

515 # Constraint-creating operators, and _predict.

516 # --------------------------------------------------------------------------

517

518 @classmethod

519 def _predict(cls, subtype, relation, other):

520 assert isinstance(subtype, cls.Subtype)

521

522 if relation == operator.__le__:

523 if issubclass(other.clstype, AffineExpression) \

524 and other.subtype.dim == 1:

525 return KullbackLeiblerConstraint.make_type(

526 argdim=subtype.argdim)

527

528 return NotImplemented

529

530 @convert_operands(scalarRHS=True)

531 @validate_prediction

532 @refine_operands()

533 def __le__(self, other):

534 if isinstance(other, AffineExpression):

535 return KullbackLeiblerConstraint(self, other)

536 else:

537 return NotImplemented

538

539

540# --------------------------------------

541__all__ = api_end(_API_START, globals())

Coverage for picos/expressions/exp_entropy.py: 61.00%

241 statements