import decimal
import os
import tempfile
from io import BytesIO
from types import MappingProxyType
import numpy as np
from getdist.paramnames import ParamInfo, ParamList, makeList
empty_dict = MappingProxyType({})
_sci_tolerance = 4
class TextFile:
def __init__(self, lines=None):
if isinstance(lines, str):
lines = [lines]
self.lines = lines or []
def write(self, outfile):
with open(outfile, "w", encoding="utf-8") as f:
f.write("\n".join(self.lines))
def texEscapeText(string):
return string.replace("_", "{\\textunderscore}")
def times_ten_power(exponent):
return r"\cdot 10^{%d}" % exponent
[docs]
def float_to_decimal(f):
# https://docs.python.org/library/decimal.html#decimal-faq
"""Convert a floating point number to a Decimal with no loss of information"""
n, d = f.as_integer_ratio()
numerator, denominator = decimal.Decimal(n), decimal.Decimal(d)
ctx = decimal.Context(prec=60)
result = ctx.divide(numerator, denominator)
while ctx.flags[decimal.Inexact]:
ctx.flags[decimal.Inexact] = False
ctx.prec *= 2
result = ctx.divide(numerator, denominator)
return result
# noinspection PyUnboundLocalVariable
def numberFigs(number, sigfig, sci=False):
# https://stackoverflow.com/questions/2663612/nicely-representing-a-floating-point-number-in-python/2663623#2663623
assert sigfig > 0
try:
d = decimal.Decimal(number)
except TypeError:
d = float_to_decimal(float(number))
if sci:
exponent = d.adjusted()
if abs(exponent) > _sci_tolerance:
d = decimal.getcontext().multiply(d, float_to_decimal(10.0**-exponent))
else:
exponent = 0
sign, digits = d.as_tuple()[0:2]
if len(digits) < sigfig:
digits = list(digits)
digits.extend([0] * (sigfig - len(digits)))
shift = d.adjusted()
result = int("".join(map(str, digits[:sigfig])))
# Round the result
if len(digits) > sigfig and digits[sigfig] >= 5:
result += 1
result = list(str(result))
# Rounding can change the length of result
# If so, adjust shift
shift += len(result) - sigfig
# reset len of result to sigfig
result = result[:sigfig]
if shift >= sigfig - 1:
# Tack more zeros on the end
result += ["0"] * (shift - sigfig + 1)
elif 0 <= shift:
# Place the decimal point in between digits
result.insert(shift + 1, ".")
else:
# Tack zeros on the front
assert shift < 0
result = ["0."] + ["0"] * (-shift - 1) + result
if sign:
result.insert(0, "-")
if sci:
return "".join(result), exponent
return "".join(result)
class NumberFormatter:
def __init__(self, sig_figs=4, separate_limit_tol=0.1, err_sf=2):
self.sig_figs = sig_figs
self.separate_limit_tol = separate_limit_tol
self.err_sf = err_sf
# noinspection PyUnboundLocalVariable
def namesigFigs(self, value, limplus, limminus, wantSign=True, sci=False):
frac = limplus / (abs(value) + limplus)
sf = self.sig_figs
if frac > 0.1 and 100 > value >= 20:
sf = 2
elif frac > 0.01 and value < 1000:
sf = 3
err_sf = self.err_sf
if value >= 20 and frac > 0.1 and limplus >= 2:
err_sf = 1
if sci:
# First, call without knowing sig figs, to get the exponent
if exponent := self.formatNumber(max(abs(value - limminus), abs(value + limplus)), sci=True)[1]:
value, limplus, limminus = (
(lambda x: decimal.getcontext().multiply(float_to_decimal(x), float_to_decimal(10.0**-exponent)))(
lim
)
for lim in [value, limplus, limminus]
)
plus_str = self.formatNumber(limplus, err_sf, wantSign)
minus_str = self.formatNumber(limminus, err_sf, wantSign)
res = self.formatNumber(value, sf)
maxdp = max(self.decimal_places(plus_str), self.decimal_places(minus_str))
while maxdp < self.decimal_places(res):
sf -= 1
if sf == 0:
res = ("%." + str(maxdp) + "f") % value
if float(res) == 0.0:
res = ("%." + str(maxdp) + "f") % 0
break
else:
res = self.formatNumber(value, sf)
while self.decimal_places(plus_str) > self.decimal_places(res):
sf += 1
res = self.formatNumber(value, sf)
if sci:
return res, plus_str, minus_str, exponent
else:
return res, plus_str, minus_str
# noinspection PyUnboundLocalVariable
def formatNumber(self, value, sig_figs=None, wantSign=False, sci=False):
if sig_figs is None:
sf = self.sig_figs
else:
sf = sig_figs
s = numberFigs(value, sf, sci=sci)
if sci:
s, exponent = s
if wantSign:
if s[0] != "-" and float(s) < 0:
s = "-" + s
if float(s) > 0:
s = "+" + s
if sci:
return s, exponent
else:
return s
def decimal_places(self, s):
i = s.find(".")
if i > 0:
return len(s) - i - 1
return 0
def plusMinusLimit(self, limit, upper, lower):
return limit != 1 or abs(abs(upper / lower) - 1) > self.separate_limit_tol
class TableFormatter:
def __init__(self):
self.border = "|"
self.endofrow = "\\\\"
self.hline = "\\hline"
self.paramText = "Parameter"
self.aboveTitles = self.hline
self.majorDividor = "|"
self.minorDividor = "|"
self.colDividor = "||"
self.belowTitles = ""
self.headerWrapper = " %s"
self.noConstraint = "---"
self.spacer = " " # just to make output more readable
self.colSeparator = self.spacer + "&" + self.spacer
self.numberFormatter = NumberFormatter()
def getLine(self, position=None):
if position is not None and hasattr(self, position):
return getattr(self, position)
return self.hline
def belowTitleLine(self, colsPerParam, numResults=None):
return self.getLine("belowTitles")
def startTable(self, ncol, colsPerResult, numResults):
part = self.majorDividor + (" c" + self.minorDividor) * (colsPerResult - 1) + " c"
return (
"\\begin{tabular} {"
+ self.border
+ " l "
+ part * numResults
+ (self.colDividor + " l " + part * numResults) * (ncol - 1)
+ self.border
+ "}"
)
def endTable(self):
return "\\end{tabular}"
def titleSubColumn(self, colsPerResult, title):
return (
" \\multicolumn{"
+ str(colsPerResult)
+ "}{"
+ self.majorDividor
+ "c"
+ self.majorDividor
+ "}{"
+ self.formatTitle(title)
+ "}"
)
def formatTitle(self, title):
return "\\bf " + texEscapeText(title)
def texEquation(self, txt):
if txt and txt[0] != "$":
return "$" + txt + "$"
else:
return txt
def textAsColumn(self, txt, latex=False, separator=False, bold=False):
wid = len(txt)
if latex:
wid += 2
if bold:
wid += 11
res = txt + self.spacer * max(0, 28 - wid)
if latex:
res = self.texEquation(res)
if bold:
res = "{\\boldmath" + res + "}"
if separator:
res += self.colSeparator
return res
class OpenTableFormatter(TableFormatter):
def __init__(self):
super().__init__()
self.border = ""
self.aboveTitles = (
r"\noalign{\vskip 3pt}" + self.hline + r"\noalign{\vskip 1.5pt}" + self.hline + r"\noalign{\vskip 5pt}"
)
self.belowTitles = r"\noalign{\vskip 3pt}" + self.hline
self.aboveHeader = ""
self.belowHeader = self.hline
self.minorDividor = ""
self.belowFinalRow = ""
def titleSubColumn(self, colsPerResult, title):
return " \\multicolumn{" + str(colsPerResult) + "}{" + "c" + "}{" + self.formatTitle(title) + "}"
class NoLineTableFormatter(OpenTableFormatter):
def __init__(self):
super().__init__()
self.aboveHeader = ""
# self.belowHeader = r'\noalign{\vskip 5pt}'
self.minorDividor = ""
self.majorDividor = ""
self.belowFinalRow = self.hline
self.belowBlockRow = self.hline
self.colDividor = "|"
self.hline = ""
def belowTitleLine(self, colsPerParam, numResults=None):
return r"\noalign{\vskip 3pt}\cline{2-" + str(colsPerParam * numResults + 1) + r"}\noalign{\vskip 3pt}"
[docs]
class ResultTable:
"""
Class for holding a latex table of parameter statistics
"""
def __init__(
self,
ncol,
results,
limit=2,
tableParamNames=None,
titles=None,
formatter=None,
numFormatter=None,
blockEndParams=None,
paramList=None,
refResults=None,
shiftSigma_indep=False,
shiftSigma_subset=False,
):
"""
:param ncol: number of columns
:param results: a :class:`MargeStats` or :class:`BestFit` instance, or a list of them for
comparing different results, or an MCSamples instance for white getMargeStats() will be called.
:param limit: which limit to include (1 is first limit calculated, usually 68%, 2 the second, usually 95%)
:param tableParamNames: optional :class:`~.paramnames.ParamNames` instance listing particular
parameters to include
:param titles: optional titles describing different results
:param formatter: a table formatting class
:param numFormatter: a number formatting class
:param blockEndParams: mark parameters in blocks, ending on this list of parameter names
:param paramList: a list of parameter names (strings) to include
:param refResults: for showing parameter shifts, a reference :class:`MargeStats` instance to show differences to
:param shiftSigma_indep: show parameter shifts in sigma assuming data are independent
:param shiftSigma_subset: show parameter shifts in sigma assuming data are a subset of each other
"""
results = list(makeList(results))
for i, res in enumerate(results):
if (getMargeStats := getattr(res, "getMargeStats", None)) is not None:
results[i] = getMargeStats()
self.lines = []
if formatter is None:
self.format = NoLineTableFormatter()
else:
self.format = formatter
self.ncol = ncol
if tableParamNames is None:
self.tableParamNames = results[0]
else:
self.tableParamNames = tableParamNames
if paramList is not None:
self.tableParamNames = self.tableParamNames.filteredCopy(paramList)
if numFormatter is not None:
self.format.numFormatter = numFormatter
self.results = results
self.boldBaseParameters = True
self.colsPerResult = len(results[0].getColumnLabels(limit))
self.colsPerParam = len(results) * self.colsPerResult
self.limit = limit
self.refResults = refResults
self.shiftSigma_indep = shiftSigma_indep
self.shiftSigma_subset = shiftSigma_subset
nparams = self.tableParamNames.numParams()
numrow = nparams // ncol
if nparams % ncol != 0:
numrow += 1
rows = []
for par in self.tableParamNames.names[0:numrow]:
rows.append([par])
for col in range(1, ncol):
for i in range(numrow * col, min(numrow * (col + 1), nparams)):
rows[i - numrow * col].append(self.tableParamNames.names[i])
self.lines.append(self.format.startTable(ncol, self.colsPerResult, len(results)))
if titles is not None:
self.addTitlesRow(titles)
self.addHeaderRow()
for row in rows[:-1]:
self.addFullTableRow(row)
if ncol == 1 and blockEndParams is not None and row[0].name in blockEndParams:
self.addLine("belowBlockRow")
else:
self.addLine("belowRow")
self.addFullTableRow(rows[-1])
self.addLine("belowFinalRow")
self.endTable()
def addFullTableRow(self, row):
txt = self.format.colSeparator.join(self.paramLabelColumn(param) + self.paramResultsTex(param) for param in row)
if not self.ncol == len(row):
txt += self.format.colSeparator * ((1 + self.colsPerParam) * (self.ncol - len(row)))
self.lines.append(txt + self.format.endofrow)
def addLine(self, position):
if self.format.getLine(position) is None: # no line is appended if the attribute is None
return self.lines
else:
return self.lines.append(self.format.getLine(position))
def addTitlesRow(self, titles):
self.addLine("aboveTitles")
cols = [self.format.titleSubColumn(1, "")]
cols += [self.format.titleSubColumn(self.colsPerResult, title) for title in titles]
self.lines.append(self.format.colSeparator.join(cols * self.ncol) + self.format.endofrow)
belowTitleLine = self.format.belowTitleLine(self.colsPerResult, self.colsPerParam // self.colsPerResult)
if belowTitleLine:
self.lines.append(belowTitleLine)
def addHeaderRow(self):
self.addLine("aboveHeader")
cols = [self.format.headerWrapper % self.format.paramText]
for result in self.results:
cols += [self.format.headerWrapper % s for s in result.getColumnLabels(self.limit)]
self.lines.append(self.format.colSeparator.join(cols * self.ncol) + self.format.endofrow)
self.addLine("belowHeader")
def paramResultsTex(self, param):
return self.format.colSeparator.join(self.paramResultTex(result, param) for result in self.results)
def paramResultTex(self, result, p):
values = result.texValues(
self.format,
p,
self.limit,
self.refResults,
shiftSigma_subset=self.shiftSigma_subset,
shiftSigma_indep=self.shiftSigma_indep,
)
if values is not None:
if len(values) > 1:
txt = self.format.textAsColumn(values[1], True, separator=True)
else:
txt = ""
txt += self.format.textAsColumn(values[0], values[0] != self.format.noConstraint)
return txt
else:
return self.format.textAsColumn("") * len(result.getColumnLabels(self.limit))
def paramLabelColumn(self, param):
return self.format.textAsColumn(param.getLabel(), True, separator=True, bold=not param.isDerived)
def endTable(self):
self.lines.append(self.format.endTable())
[docs]
def tableTex(self, document=False, latex_preamble=None, packages=("amsmath", "amssymb", "bm")):
"""
Get the latex string for the table
:param document: if True, make a full latex file, if False just the snippet for including in another file
:param latex_preamble: any preamble to include in the latex file
:param packages: list of packages to load
"""
if document:
lines = [r"\documentclass{article}", r"\pagestyle{empty}"]
for package in packages:
lines.append(r"\usepackage{%s}" % package)
lines.append("\\renewcommand{\\arraystretch}{1.5}")
if latex_preamble:
lines.append(latex_preamble)
lines.append("\\begin{document}")
lines += self.lines
lines.append("\\end{document}")
else:
lines = self.lines
return "\n".join(lines)
[docs]
def write(self, fname, **kwargs):
"""
Write the latex for the table to a file
:param fname: filename to write
:param kwargs: arguments for :func:`~ResultTable.tableTex`
"""
TextFile(self.tableTex(**kwargs)).write(fname)
[docs]
def tablePNG(self, dpi=None, latex_preamble=None, filename=None, bytesIO=False):
"""
Get a .png file image of the table. You must have latex installed to use this.
:param dpi: dpi settings for the png
:param latex_preamble: any latex preamble
:param filename: filename to save to (defaults to file in the temp directory)
:param bytesIO: if True, return a BytesIO instance holding the .png data
:return: if bytesIO, the BytesIO instance, otherwise name of the output file
"""
import subprocess
texfile = tempfile.mktemp(suffix=".tex")
self.write(texfile, document=True, latex_preamble=latex_preamble)
basefile = os.path.splitext(texfile)[0]
outfile = filename or basefile + ".png"
old_pwd = os.getcwd()
def runCommand(command):
# Use CREATE_NO_WINDOW flag on Windows to prevent console window from appearing
creationflags = 0
if os.name == "nt": # Windows
creationflags = subprocess.CREATE_NO_WINDOW
try:
subprocess.run(
command,
stdout=subprocess.DEVNULL,
stderr=subprocess.DEVNULL,
creationflags=creationflags,
check=True, # Raise exception on non-zero return code
)
except FileNotFoundError:
cmd_name = command[0] if command else "Command"
error_msg = f"Command not found: {cmd_name}"
if cmd_name == "latex":
error_msg += "\nLaTeX must be installed to generate tables. Please install a TeX distribution like TeX Live, MiKTeX, or MacTeX."
elif cmd_name == "dvipng":
error_msg += (
"\ndvipng must be installed to generate PNG images. It is included in most LaTeX distributions."
)
raise FileNotFoundError(error_msg)
try:
os.chdir(os.path.dirname(texfile))
runCommand(["latex", texfile])
cmd = ["dvipng"]
if dpi:
cmd.extend(["-D", str(dpi)])
cmd.extend(["-T", "tight", "-x", "1000", "-z", "9", "--truecolor", "-o", outfile, basefile + ".dvi"])
runCommand(cmd)
finally:
for f in [basefile + ext for ext in (".tex", ".dvi", ".aux", ".log")]:
if os.path.isfile(f):
os.remove(f)
os.chdir(old_pwd)
if bytesIO:
with open(outfile, "rb") as f:
result = BytesIO(f.read())
os.remove(outfile)
result.seek(0)
return result
else:
return outfile
[docs]
class ParamResults(ParamList):
"""
Base class for a set of parameter results, inheriting from :class:`~.paramnames.ParamList`,
so that self.names is a list of :class:`~.paramnames.ParamInfo` instances for each parameter, which
have attribute holding results for the different parameters.
"""
class LikelihoodChi2:
name: str
tag: str
chisq: float
[docs]
class BestFit(ParamResults):
"""
Class holding the result of a likelihood minimization, inheriting from :class:`ParamResults`.
The data is read from a specific formatted text file (.minimum or .bestfit) as output by CosmoMC or Cobaya.
"""
def __init__(self, fileName=None, setParamNameFile=None, want_fixed=False, max_posterior=True):
"""
:param fileName: text file to load from, assumed to be in CosmoMC's .minimum format
:param setParamNameFile: optional name of .paramnames file listing preferred parameter labels for the parameters
:param want_fixed: whether to include values of parameters that are not allowed to vary
:param max_posterior: whether the file is a maximum posterior (default) or maximum likelihood
"""
super().__init__()
self.max_posterior = max_posterior
if fileName is not None:
self.loadFromFile(fileName, want_fixed=want_fixed)
if setParamNameFile is not None:
self.setLabelsFromParamNames(setParamNameFile)
def getColumnLabels(self, **_kwargs):
return ["Best fit"]
def loadFromFile(self, filename, want_fixed=False):
textFileLines = self.fileList(filename)
first = textFileLines[0].strip().split("=")
if first[0].strip() == "weight":
self.weight = float(first[1].strip())
del textFileLines[0]
first = textFileLines[0].strip().split("=")
if first[0].strip() != "-log(Like)":
raise Exception("Error in format of parameter (best fit) file")
self.logLike = float(first[1].strip())
isFixed = False
isDerived = False
self.chiSquareds = []
chunks = 0
if len(textFileLines[1].strip()) > 0:
del textFileLines[1] # if it has chi2 line as well
for ix in range(2, len(textFileLines)):
line = textFileLines[ix]
if len(line.strip()) == 0:
chunks += 1
isFixed = not isFixed
isDerived = True
if chunks == 3:
if ix + 2 >= len(textFileLines):
break
for likePart in textFileLines[ix + 2 :]:
if len(likePart.strip()) != 0:
(chisq, name) = [s.strip() for s in likePart.split(None, 2)][1:]
name = [s.strip() for s in name.split(":", 1)]
if len(name) > 1:
(kind, name) = name
else:
kind = ""
chi2 = LikelihoodChi2()
if "=" in name:
chi2.tag, chi2.name = (s.strip() for s in name.split("="))
else:
chi2.tag, chi2.name = None, name
chi2.chisq = float(chisq)
self.chiSquareds.append((kind, chi2))
break
continue
if not isFixed or want_fixed:
param = ParamInfo()
param.isFixed = isFixed
param.isDerived = isDerived
(param.number, param.best_fit, param.name, param.label) = (s.strip() for s in line.split(None, 3))
param.number = int(param.number)
param.best_fit = float(param.best_fit)
self.names.append(param)
def sortedChiSquareds(self):
likes = dict()
for kind, val in self.chiSquareds:
if kind not in likes:
likes[kind] = []
likes[kind].append(val)
return sorted(iter(likes.items()))
def chiSquareForKindName(self, kind, name):
for akind, val in self.chiSquareds:
if akind == kind and val.name == name:
return val.chisq
return None
def texValues(self, formatter, p, **_kwargs):
param = self.parWithName(p.name)
if param is not None:
return [formatter.numberFormatter.formatNumber(param.best_fit)]
else:
return None
def getParamDict(self, include_derived=True):
res = dict()
for i, name in enumerate(self.names):
if include_derived or not name.isDerived:
res[name.name] = name.best_fit
res["weight"] = 1
res["loglike"] = self.logLike
return res
[docs]
class ParamLimit:
"""
Class containing information about a marginalized parameter limit.
:ivar lower: lower limit
:ivar upper: upper limit
:ivar twotail: True if a two-tail limit, False if one-tail
:ivar onetail_upper: True if one-tail upper limit
:ivar onetail_lower: True if one-tail lower limit
"""
def __init__(self, minmax, tag="two"):
"""
:param minmax: a [min,max] tuple with lower and upper limits. Entries be None if no limit.
:param tag: a text tag descibing the limit, one of ['two' | '>' | '<' | 'none']
"""
self.lower = minmax[0]
self.upper = minmax[1]
self.twotail = tag == "two"
self.onetail_upper = tag == ">"
self.onetail_lower = tag == "<"
[docs]
def limitTag(self):
"""
:return: Short text tag describing the type of limit (one-tail or two tail):
- *two*: two-tail limit
- *>*: a one-tail upper limit
- *<*: a one-tail lower limit
- *none*: no limits (both boundaries have high probability)
"""
if self.twotail:
return "two"
elif self.onetail_upper:
return ">"
elif self.onetail_lower:
return "<"
else:
return "none"
[docs]
def limitType(self):
"""
:return: a text description of the type of limit. One of:
- *two tail*
- *one tail upper limit*
- *one tail lower limit*
- *none*
"""
if self.twotail:
return "two tail"
elif self.onetail_upper:
return "one tail upper limit"
elif self.onetail_lower:
return "one tail lower limit"
else:
return "none"
def __str__(self):
"""
:return: string representation of lower and upper bounds, with text description of the limit type
"""
return f"{self.lower:g} {self.upper:g} {self.limitTag()}"
[docs]
class MargeStats(ParamResults):
"""
Stores marginalized 1D parameter statistics, including mean, variance and confidence limits,
inheriting from :class:`ParamResults`.
Values are stored as attributes of the :class:`~.paramnames.ParamInfo` objects stored in self.names.
Use *par= margeStats.parWithName('xxx')* to get the :class:`~.paramnames.ParamInfo` for parameter *xxx*;
Values stored are:
- *par.mean*: parameter mean
- *par.err*: standard deviation
- *limits*: list of :class:`~.types.ParamLimit` objects for the stored number of marginalized limits
For example to get the first and second lower limits (default 68% and 95%) for parameter *xxx*::
print(margeStats.names.parWithName('xxx').limits[0].lower)
print(margeStats.names.parWithName('xxx').limits[1].lower)
See :class:`~.types.ParamLimit` for details of limits.
"""
[docs]
def loadFromFile(self, filename):
"""
Load from a plain text file
:param filename: file to load from
"""
textFileLines = self.fileList(filename)
lims = textFileLines[0].split(":")[1]
self.limits = [float(s.strip()) for s in lims.split(";")]
self.hasBestFit = False
for line in textFileLines[3:]:
if len(line.strip()) == 0:
break
param = ParamInfo()
items = [s.strip() for s in line.split(None, len(self.limits) * 3 + 3)]
param.name = items[0]
if param.name[-1] == "*":
param.isDerived = True
param.name = param.name[:-1]
param.mean = float(items[1])
param.err = float(items[2])
param.label = items[-1]
param.limits = []
for i in range(len(self.limits)):
param.limits.append(ParamLimit([float(s) for s in items[3 + i * 3 : 5 + i * 3]], items[5 + i * 3]))
self.names.append(param)
def headerLine(self, inc_limits=False):
parForm = self.parFormat()
text = ""
text += parForm % "parameter" + " "
text += "%-15s" % "mean"
text += "%-15s" % "sddev"
for j, limit in enumerate(self.limits):
if inc_limits:
tag = "_%.0f%%" % (limit * 100)
limtxt = "type"
else:
tag = str(j + 1)
limtxt = "limit" + tag
text += "%-15s" % ("lower" + tag)
text += "%-15s" % ("upper" + tag)
text += "%-7s" % limtxt
return text, parForm
def __str__(self):
contours_str = "; ".join([str(c) for c in self.limits])
header, parForm = self.headerLine()
text = ""
text += "Marginalized limits: %s\n\n" % contours_str
text += header
text += "\n"
for j, par in enumerate(self.names):
text += parForm % (self.name(j, True))
text += f"{par.mean:15.7E}{par.err:15.7E}"
for lim in par.limits:
text += "%15.7E%15.7E %-5s" % (lim.lower, lim.upper, lim.limitTag())
text += " %s\n" % par.label
return text
def addBestFit(self, bf):
self.hasBestFit = True
self.logLike = bf.logLike
# the next line deletes parameters not in best-fit;
# this is good e.g. to get rid of yhe from importance sampled result
self.names = [x for x in self.names if bf.parWithName(x.name) is not None]
for par in self.names:
param = bf.parWithName(par.name)
par.best_fit = param.best_fit
par.isDerived = param.isDerived
def limitText(self, limit):
txt = str(round(self.limits[limit - 1] * 100.0))
if txt.endswith(".0"): # e.g. 95.0 -> 95
txt = txt.split(".")[0]
return txt
def getColumnLabels(self, limit=2):
if self.hasBestFit:
res = ["Best fit"]
else:
res = []
return res + [self.limitText(limit) + "\\% limits"]
def texValues(self, formatter, p, limit=2, refResults=None, shiftSigma_indep=False, shiftSigma_subset=False):
if not isinstance(p, ParamInfo):
param = self.parWithName(p)
else:
param = self.parWithName(p.name)
if param is not None:
lim = param.limits[limit - 1]
sf = 3
if param.name.startswith("chi2"):
# Chi2 for low dof are very skewed, always want mean and sigma or limit
res, sigma, _ = formatter.numberFormatter.namesigFigs(
param.mean, param.err, param.err, wantSign=False, sci=False
)
if limit == 1:
res += r"\pm " + sigma
else:
# in this case give mean and effective dof
res += r"\,({\nu\rm{:}\,%.1f})" % (param.err**2 / 2)
# res, plus_str, minus_str =
# formatter.numberFormatter.namesigFigs(param.mean, lim.upper - param.mean, lim.lower, sci=False)
# res += '^{' + plus_str + '}_{>' + minus_str + '}'
elif lim.twotail:
if not formatter.numberFormatter.plusMinusLimit(limit, lim.upper - param.mean, lim.lower - param.mean):
res, plus_str, _, exponent = formatter.numberFormatter.namesigFigs(
param.mean, param.err, param.err, wantSign=False, sci=True
)
res += r"\pm " + plus_str
else:
res, plus_str, minus_str, exponent = formatter.numberFormatter.namesigFigs(
param.mean, lim.upper - param.mean, lim.lower - param.mean, sci=True
)
res += "^{" + plus_str + "}_{" + minus_str + "}"
if exponent:
res = r"\left(\,%s\,\right)" % res + times_ten_power(exponent)
elif lim.onetail_upper:
res, exponent = formatter.numberFormatter.formatNumber(lim.upper, sf, sci=True)
res = "< " + res
if exponent:
res += times_ten_power(exponent)
elif lim.onetail_lower:
res, exponent = formatter.numberFormatter.formatNumber(lim.lower, sf, sci=True)
res = "> " + res
if exponent:
res += times_ten_power(exponent)
else:
res = formatter.noConstraint
if refResults is not None and res != formatter.noConstraint:
refVal = refResults.parWithName(param.name)
if refVal is not None:
delta = param.mean - refVal.mean
if shiftSigma_indep or shiftSigma_subset:
res += r"\quad("
if shiftSigma_subset:
# give mean shift in sigma units for subset data (regularized to max sigma/20)
subset_sigma = np.sqrt(abs(param.err**2 - refVal.err**2))
res += "%+.1f \\sigma_s" % (delta / max(subset_sigma, refVal.err / 20))
if shiftSigma_indep:
# give mean shift in sigma units for independent data
indep_sigma = np.sqrt(param.err**2 + refVal.err**2)
res += ", %+.1f \\sigma_i" % (delta / indep_sigma)
res += ")"
else:
res += r"\quad(%+.1f \sigma)" % (delta / refVal.err)
if self.hasBestFit: # add best fit too
rangew = (lim.upper - lim.lower) / 10
bestfit, _, _, exponent = formatter.numberFormatter.namesigFigs(
param.best_fit, rangew, -rangew, sci=True
)
if exponent:
bestfit += times_ten_power(exponent)
return [res, bestfit]
return [res]
else:
return None
[docs]
class LikeStats(ParamResults):
"""
Stores posterior-related statistics, including best-fit sample and extremal values of the N-D confidence region,
inheriting from :class:`ParamResults`. Note: despite the name "LikeStats", this actually stores statistics
related to the posterior (likelihood × prior), not just the likelihood.
TODO: currently only saves to text, does not load full data from file
"""
def loadFromFile(self, filename):
textFileLines = self.fileList(filename)
results = dict()
for line in textFileLines:
if len(line.strip()) == 0:
break
name, value = (x.strip() for x in line.split("="))
results[name] = float(value)
self.logLike_sample = results.get("Best fit sample -log(Like)")
self.logMeanInvLike = results.get("Ln(mean 1/like)")
self.meanLogLike = results.get("mean(-Ln(like))")
self.logMeanLike = results.get("-Ln(mean like)")
self.complexity = results.get("complexity")
twiceVarLogLike = results.get("2*Var(Ln(like))")
if twiceVarLogLike is not None:
self.varLogLike = 0.5 * twiceVarLogLike
else:
self.varLogLike = None
# TODO: load N-D limits
def likeSummary(self):
text = "Best fit sample -log(Like) = %f\n" % self.logLike_sample
if self.logMeanInvLike:
text += "Ln(mean 1/like) = %f\n" % self.logMeanInvLike
text += "mean(-Ln(like)) = %f\n" % self.meanLogLike
text += "-Ln(mean like) = %f\n" % self.logMeanLike
text += "2*Var(Ln(like)) = %f\n" % (self.varLogLike * 2.0)
# text += "complexity = %f\n" % self.complexity
return text
def headerLine(self):
return self.parFormat() % "parameter" + " bestfit lower1 upper1 lower2 upper2\n"
def __str__(self):
text = self.likeSummary()
parForm = self.parFormat()
if len(self.names):
text += "\n"
text += self.headerLine()
for j, par in enumerate(self.names):
if par.ND_limit_bot.size < 2:
raise Exception("Likestats output assumes at least two contour levels")
text += parForm % (self.name(j, True))
text += f"{par.bestfit_sample:15.7E}{par.ND_limit_bot[0]:15.7E}{par.ND_limit_top[0]:15.7E}{par.ND_limit_bot[1]:15.7E}{par.ND_limit_top[1]:15.7E} {par.label}\n"
return text
[docs]
class ConvergeStats(ParamResults):
def loadFromFile(self, filename):
try:
textFileLines = self.fileList(filename)
self.R_eigs = []
for i in range(len(textFileLines)):
if textFileLines[i].find("var(mean)") >= 0:
for line in textFileLines[i + 1 :]:
if len(line.strip()) == 0:
break
try:
self.R_eigs.append(line.split()[1])
except Exception:
self.R_eigs.append("1e30")
elif "Parameter auto-correlations" in textFileLines[i]:
self.auto_correlation_steps = [int(s) for s in textFileLines[i + 2].split()]
self.auto_correlations = []
self.auto_correlation_pars = []
for line in textFileLines[i + 3 :]:
if len(line.strip()) == 0:
break
items = line.split(None, len(self.auto_correlation_steps) + 1)
self.auto_correlation_pars.append(items[0])
self.auto_correlations.append([float(s) for s in items[1:-1]])
except:
print("Error reading: " + filename)
raise
def worstR(self, default=None):
return self.R_eigs[len(self.R_eigs) - 1] if self.R_eigs else default