From 2df8bf24439610249edd01b99eb281a7f85e2d17 Mon Sep 17 00:00:00 2001 From: beguinotj Date: Mon, 6 Oct 2025 14:22:34 +0200 Subject: [PATCH 01/12] Code, test and doc for MGL --- docs/refs.bib | 22 ++ src/scalib/postprocessing/__init__.py | 4 +- .../postprocessing/noise_amplification.py | 188 ++++++++++++++++++ 3 files changed, 213 insertions(+), 1 deletion(-) create mode 100644 src/scalib/postprocessing/noise_amplification.py diff --git a/docs/refs.bib b/docs/refs.bib index 5fbb0f0b..fda4f173 100644 --- a/docs/refs.bib +++ b/docs/refs.bib @@ -89,3 +89,25 @@ @inproceedings{rankest_histograms publisher = {Springer}, year = {2016} } + +@inproceedings{mgl, + author = {Julien B{\'{e}}guinot and + Wei Cheng and + Sylvain Guilley and + Yi Liu and + Lo{\"{\i}}c Masure and + Olivier Rioul and + Fran{\c{c}}ois{-}Xavier Standaert}, + editor = {Elif Bilge Kavun and + Michael Pehl}, + title = {Removing the Field Size Loss from Duc et al.'s Conjectured Bound for + Masked Encodings}, + booktitle = {Constructive Side-Channel Analysis and Secure Design - 14th International + Workshop, {COSADE} 2023, Munich, Germany, April 3-4, 2023, Proceedings}, + series = {Lecture Notes in Computer Science}, + volume = {13979}, + pages = {86--104}, + publisher = {Springer}, + year = {2023} +} + diff --git a/src/scalib/postprocessing/__init__.py b/src/scalib/postprocessing/__init__.py index 85b2cf0d..4f40e4af 100644 --- a/src/scalib/postprocessing/__init__.py +++ b/src/scalib/postprocessing/__init__.py @@ -10,8 +10,10 @@ :nosignatures: scalib.postprocessing.rankestimation + scalib.postprocessing.noise_amplification """ -__all__ = ["rankestimation"] +__all__ = ["rankestimation", "mrs_gerber_lemma"] from .rankestimation import rank_nbin, rank_accuracy +from .noise_amplification import mrs_gerber_lemma diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py new file mode 100644 index 00000000..dfaeeeb1 --- /dev/null +++ b/src/scalib/postprocessing/noise_amplification.py @@ -0,0 +1,188 @@ +r"""Estimation of the mutual information between a sensitive value protected by masking and leakages in terms of the mutual information between each share and its corresponding leakages + +The `mrs_gerber_lemma` function takes as input the mutual information for each share separately (eventually for multiple sensitive values) +and outputs an upper upper bound on the mutual information between the sensitive value and the leakages. + +By default, it is assumed that the leakages are expressed in bits. +Eventually, a specific base for the unit of information can be specified. + +The derivation is based on Mrs Gerber's lemma. +In particular, it assumes that the the shares are leaking separately which should be ensured by a proper implementation of the masking countermeasure. + +Examples +-------- + +>>> from scalib.postprocessing import mrs_gerber_lemma +>>> import numpy as np +>>> # Mutual information on three shares +>>> MI_shares = np.array([.1,.2,.5]) +>>> # Derive an upper bound on the mutual information of the protected secret +>>> MI_sensitive = mrs_gerber_lemma(MI_shares,group_order=2**8, base=2) + +Reference +--------- + +.. currentmodule:: scalib.postprocessing.noise_amplification + +.. autosummary:: + :toctree: + :nosignatures: + :recursive: + + mrs_gerber_lemma + +Notes +----- +The upper bound is based on the article :footcite:p:`mgl`, + +References +---------- + +.. footbibliography:: +""" + +__all__ = ["mrs_gerber_lemma"] +import numpy as np + + +def phi(x: float) -> float: + r"""Compute the DFT of the binary entropy. + + Parameters + ---------- + x : array_like, f64 + Array of floats assumed to belong to [0,1] + + Returns + ------- + An array corresponding to the image of the input array by the DFT of the binary entropy (in nats) + """ + x = np.asarray(x) + + # Deal with special case 'x=1' with where + y = np.where( + x == 1, np.log(2), ((1 - x) * np.log1p(-x) + (1 + x) * np.log1p(x)) / 2 + ) + return y + + +def phi_derivative(x: float) -> float: + r"""Compute the derivative of the DFT of the binary entropy. + + Parameters + ---------- + x : array_like, f64 + Array of floats assumed to belong to [0,1] + + Returns + ------- + An array corresponding to the image of the input array by the derivative of the DFT of the binary entropy + """ + x = np.asarray(x) + return np.arctanh(x) + + +def phi_second_derivative(x: float) -> float: + r"""Compute the second derivative of the DFT of the binary entropy. + + Parameters + ---------- + x : array_like, f64 + Array of floats assumed to belong to [0,1] + + Returns + ------- + An array corresponding to the image of the input array by the second derivative of the DFT of the binary entropy + """ + x = np.asarray(x) + return (1 - x**2) ** (-1) + + +def phi_inv(y: float, niter=3) -> float: + r"""Compute the inverse of the DFT of the binary entropy via Halley's root finding algorithm + + Parameters + ---------- + y : array_like, f64 + Array of floats to apply the inverse, the values should belong to the interval [0,\log 2] in nats + + niter: int + Number of iteration used in Halley's method. By default niter=3. + + Returns + ------- + An array corresponding to the image of the input array by the inverse of the DFT of the binary entropy. + """ + + y = np.asarray(y) + + # Initial Guess of the inverse + # Based on https://math.stackexchange.com/questions/3454390/find-the-approximation-of-the-inverse-of-binary-entropy-function + x = np.sqrt(1 - (1 - y / np.log(2)) ** (4 / 3)) + + # Iterative Halley's method for root finding with cubic convergence rate + # See https://en.wikipedia.org/wiki/Halley's_method + for _ in range(niter): + + f = phi(x) - y + df = phi_derivative(x) + ddf = phi_second_derivative(x) + + x -= (f * df) / (df**2 - f * ddf / 2) + + # We deal with the special case 0 and log(2) with where + return np.where(y > 0, np.where(y < np.log(2), x, 1), 0) + + +def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: + r"""Upper bound the mutual information of a sensitive value in terms of the mutual information for each of its share. + + Parameters + ---------- + MI_shares : array_like, f64 + Mutual information for each share. Array must be of shape ``(ns,nv)`` where + ``ns`` is the number of shares, ``nv`` the number of sensitive values. + group_order : int + Order of the Abelian in which the sensitive values are protected by masking. + base : array_like, f64 + The base of information used, by default the information is in bits i.e. base=2. + + Returns + ------- + Upper bound on the mutual information for all nv sensitive values based on Mrs Gerber's Lemma + """ + + MI_shares = np.asarray(MI_shares) + + # Check if the group order is a power of 2 (i.e. its bit expression contains a single 1) + if (group_order & (group_order - 1)) == 0: + MI_share_bits = MI_shares * np.log2(base) + k = np.min(np.floor(MI_share_bits), axis=0) + cliped_MI_shares = np.clip(0, 1, MI_share_bits - k) * np.log(2) + MI = k * np.log(2) + phi(np.prod(phi_inv(cliped_MI_shares), axis=0)) + # Otherwise, use a weaker MGL based on Pinsker/reverse Pinsker inequalities + else: + beta = group_order**2 * 4 ** (1 / group_order) + + # Detect which shares are bellow the noise amplification ratio + bellow = 2 * MI_shares * np.log(base) < 1 + + # Multipy only the terms bellow the noise amplification ratio + product = np.prod(np.where(bellow, 2 * MI_shares * np.log(base), 1), axis=0) / 4 + + # If there is no share bellow the noise amplification ratio we return the minimum instead + min_shares = np.min(MI_shares * np.log(base), axis=0) + + # Depending on the case we return either the minimum or the product of shares bellow noise amplification + P = np.where(~np.any(bellow, axis=0), min_shares, product) + + MI_1 = np.log(1 + beta * P) + MI_2 = (1 / group_order + np.sqrt(P)) * np.log(1 + group_order * np.sqrt(P)) + min_MI = np.minimum(np.minimum(MI_1, MI_2), min_shares) + + # Depending on the case we return either the minimum or the amplification lemma + MI = np.where(~np.any(bellow, axis=0), min_shares, min_MI) + + # Conversion from nats to base 'base' + MI /= np.log(base) + return MI From 4a01043d93961dbf7562ecf60af8b9fe07a9ae1b Mon Sep 17 00:00:00 2001 From: beguinotj Date: Mon, 6 Oct 2025 14:33:38 +0200 Subject: [PATCH 02/12] Missing test file --- tests/test_mgl.py | 33 +++++++++++++++++++++++++++++++++ 1 file changed, 33 insertions(+) create mode 100644 tests/test_mgl.py diff --git a/tests/test_mgl.py b/tests/test_mgl.py new file mode 100644 index 00000000..a4bef3bc --- /dev/null +++ b/tests/test_mgl.py @@ -0,0 +1,33 @@ +import pytest +from scalib.postprocessing import mrs_gerber_lemma +import numpy as np + + +def test_mgl(): + + err = 10**-15 + x = np.random.uniform(low=0, high=1) + + # Test with a single share + assert abs(mrs_gerber_lemma(x, 2**8, base=2) - x) <= err + + # Test when a share leak more than a bit and another one less than a bit + assert abs(mrs_gerber_lemma([x, 1.1], 2**8, base=2) - x) <= err + + # Test when two shares leaks between 1 and 2 bits + assert ( + abs( + mrs_gerber_lemma([1 + x, 1.1], 2**8, base=2) + - (1 + mrs_gerber_lemma([x, 0.1], 2**8, base=2)) + ) + <= err + ) + + # Test when a single share is bellow noise amplification + assert abs(mrs_gerber_lemma(x, 2**8 - 1, base=2) - x) <= err + + # Test when one share is bellow noise amplification threshold + assert abs(mrs_gerber_lemma([x, 15], 2**8 - 1, base=2) - x) <= err + + # Test when no share is bellow noise amplification threshold + assert abs(mrs_gerber_lemma([12, 15], 2**8 - 1, base=2) - 12) <= err From 7370c873d93ca20774c9bdc14c632bc68c110f0b Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:22:44 +0200 Subject: [PATCH 03/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index dfaeeeb1..1a25cf8b 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -143,7 +143,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: Mutual information for each share. Array must be of shape ``(ns,nv)`` where ``ns`` is the number of shares, ``nv`` the number of sensitive values. group_order : int - Order of the Abelian in which the sensitive values are protected by masking. + Order of the group in which the sensitive values are protected by masking. base : array_like, f64 The base of information used, by default the information is in bits i.e. base=2. From 3724db5e50f5600e5c401474cbd5523b1bb00af2 Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:22:54 +0200 Subject: [PATCH 04/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 1a25cf8b..cf7e4492 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -149,7 +149,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: Returns ------- - Upper bound on the mutual information for all nv sensitive values based on Mrs Gerber's Lemma + Upper bound on the mutual information for all nv sensitive values based on Mrs Gerber's Lemma. """ MI_shares = np.asarray(MI_shares) From 64ea6b6b4b6a0c8d262ba7b31922cbdb92487906 Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:23:08 +0200 Subject: [PATCH 05/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index cf7e4492..6d20deac 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -164,7 +164,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: else: beta = group_order**2 * 4 ** (1 / group_order) - # Detect which shares are bellow the noise amplification ratio + # Detect which shares are below the noise amplification ratio. bellow = 2 * MI_shares * np.log(base) < 1 # Multipy only the terms bellow the noise amplification ratio From 95956f5f07c6fa01b143d288217d2318749bdfae Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:23:18 +0200 Subject: [PATCH 06/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 6d20deac..207c4a8e 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -168,7 +168,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: bellow = 2 * MI_shares * np.log(base) < 1 # Multipy only the terms bellow the noise amplification ratio - product = np.prod(np.where(bellow, 2 * MI_shares * np.log(base), 1), axis=0) / 4 + product = np.prod(np.where(below, 2 * MI_shares * np.log(base), 1), axis=0) / 4 # If there is no share bellow the noise amplification ratio we return the minimum instead min_shares = np.min(MI_shares * np.log(base), axis=0) From d441de8541ede1d5087bd00c5d9ce6202c36578e Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:23:25 +0200 Subject: [PATCH 07/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 207c4a8e..63a9d7e5 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -181,7 +181,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: min_MI = np.minimum(np.minimum(MI_1, MI_2), min_shares) # Depending on the case we return either the minimum or the amplification lemma - MI = np.where(~np.any(bellow, axis=0), min_shares, min_MI) + MI = np.where(~np.any(below, axis=0), min_shares, min_MI) # Conversion from nats to base 'base' MI /= np.log(base) From 24ec43c366bb7d01c62fac6e21813be92a0625bf Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:23:33 +0200 Subject: [PATCH 08/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 63a9d7e5..23e02db9 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -173,8 +173,8 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: # If there is no share bellow the noise amplification ratio we return the minimum instead min_shares = np.min(MI_shares * np.log(base), axis=0) - # Depending on the case we return either the minimum or the product of shares bellow noise amplification - P = np.where(~np.any(bellow, axis=0), min_shares, product) + # Depending on the case we return either the minimum or the product of shares below noise amplification + P = np.where(~np.any(below, axis=0), min_shares, product) MI_1 = np.log(1 + beta * P) MI_2 = (1 / group_order + np.sqrt(P)) * np.log(1 + group_order * np.sqrt(P)) From a0ad5e7d13c4a59a9c424873bffb1ba3b53d915f Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:23:45 +0200 Subject: [PATCH 09/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 23e02db9..aff62639 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -158,8 +158,8 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: if (group_order & (group_order - 1)) == 0: MI_share_bits = MI_shares * np.log2(base) k = np.min(np.floor(MI_share_bits), axis=0) - cliped_MI_shares = np.clip(0, 1, MI_share_bits - k) * np.log(2) - MI = k * np.log(2) + phi(np.prod(phi_inv(cliped_MI_shares), axis=0)) + clipped_MI_shares = np.clip(0, 1, MI_share_bits - k) * np.log(2) + MI = k * np.log(2) + phi(np.prod(phi_inv(clipped_MI_shares), axis=0)) # Otherwise, use a weaker MGL based on Pinsker/reverse Pinsker inequalities else: beta = group_order**2 * 4 ** (1 / group_order) From 7ad05a70277e49372bd2ac72d66be2cb41541e5c Mon Sep 17 00:00:00 2001 From: JulienBeg <57485155+JulienBeg@users.noreply.github.com> Date: Wed, 15 Oct 2025 23:24:05 +0200 Subject: [PATCH 10/12] Update src/scalib/postprocessing/noise_amplification.py MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Gaëtan Cassiers --- src/scalib/postprocessing/noise_amplification.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index aff62639..aca57595 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -165,7 +165,7 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: beta = group_order**2 * 4 ** (1 / group_order) # Detect which shares are below the noise amplification ratio. - bellow = 2 * MI_shares * np.log(base) < 1 + below = 2 * MI_shares * np.log(base) < 1 # Multipy only the terms bellow the noise amplification ratio product = np.prod(np.where(below, 2 * MI_shares * np.log(base), 1), axis=0) / 4 From 6206a857f94651923ce276879f3609ab498b53db Mon Sep 17 00:00:00 2001 From: beguinotj Date: Sun, 19 Oct 2025 18:06:06 +0200 Subject: [PATCH 11/12] Modifications based on the code review --- src/scalib/postprocessing/__init__.py | 4 +- .../postprocessing/noise_amplification.py | 105 ++++++++++++------ tests/test_mgl.py | 45 ++++++-- 3 files changed, 103 insertions(+), 51 deletions(-) diff --git a/src/scalib/postprocessing/__init__.py b/src/scalib/postprocessing/__init__.py index 4f40e4af..bb144a19 100644 --- a/src/scalib/postprocessing/__init__.py +++ b/src/scalib/postprocessing/__init__.py @@ -13,7 +13,7 @@ scalib.postprocessing.noise_amplification """ -__all__ = ["rankestimation", "mrs_gerber_lemma"] +__all__ = ["rankestimation", "mgl"] from .rankestimation import rank_nbin, rank_accuracy -from .noise_amplification import mrs_gerber_lemma +from .noise_amplification import mgl diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index aca57595..14a02066 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -1,6 +1,13 @@ r"""Estimation of the mutual information between a sensitive value protected by masking and leakages in terms of the mutual information between each share and its corresponding leakages -The `mrs_gerber_lemma` function takes as input the mutual information for each share separately (eventually for multiple sensitive values) +This function is usefull in the following setting. +You know that a sensitive value X is protected by masking so that it is shared into (S_0,...,S_d). +You can observe leakages Y_0,...,Y_d for each corresponding shares. +You have have estimated the leakages on each share via the mutual information I(S_i;Y_i). +Then, the mgl function provides an upper bound on the mutual information I(X; Y_0,...,Y_d). +The obtained upper bound can then be used with other functions that provides security guarantees (such as success rate of an attack) in terms of mutual information. + +The `mgl` function takes as input the mutual information for each share separately (eventually for multiple sensitive values) and outputs an upper upper bound on the mutual information between the sensitive value and the leakages. By default, it is assumed that the leakages are expressed in bits. @@ -12,12 +19,12 @@ Examples -------- ->>> from scalib.postprocessing import mrs_gerber_lemma +>>> from scalib.postprocessing import mgl >>> import numpy as np >>> # Mutual information on three shares ->>> MI_shares = np.array([.1,.2,.5]) +>>> mi_shares = np.array([.1,.2,.5]) >>> # Derive an upper bound on the mutual information of the protected secret ->>> MI_sensitive = mrs_gerber_lemma(MI_shares,group_order=2**8, base=2) +>>> mi_sensitive = mgl(mi_shares,group_order=2**8, base=2) Reference --------- @@ -29,7 +36,7 @@ :nosignatures: :recursive: - mrs_gerber_lemma + mgl Notes ----- @@ -41,13 +48,19 @@ .. footbibliography:: """ -__all__ = ["mrs_gerber_lemma"] +__all__ = ["mgl"] import numpy as np +import numpy.typing as npt -def phi(x: float) -> float: +def phi(x: npt.NDArray[np.float64]) -> npt.NDArray[np.float64]: r"""Compute the DFT of the binary entropy. + See equation 51 (Theorem 1) in the book chapter: + Olivier Rioul, Julien Béguinot. The role of Mrs. Gerber’s Lemma for evaluating the information + leakage of secret sharing schemes. Ioannis Kontoyiannis, Jason Klusowski, Cynthia Rush. Information + Theory, Probability and Statistical Learning: A Festschrift in Honor of Andrew Barron, Springer, 2025. + Parameters ---------- x : array_like, f64 @@ -57,8 +70,6 @@ def phi(x: float) -> float: ------- An array corresponding to the image of the input array by the DFT of the binary entropy (in nats) """ - x = np.asarray(x) - # Deal with special case 'x=1' with where y = np.where( x == 1, np.log(2), ((1 - x) * np.log1p(-x) + (1 + x) * np.log1p(x)) / 2 @@ -66,7 +77,7 @@ def phi(x: float) -> float: return y -def phi_derivative(x: float) -> float: +def phi_derivative(x: npt.NDArray[np.float64]) -> npt.NDArray[np.float64]: r"""Compute the derivative of the DFT of the binary entropy. Parameters @@ -78,11 +89,10 @@ def phi_derivative(x: float) -> float: ------- An array corresponding to the image of the input array by the derivative of the DFT of the binary entropy """ - x = np.asarray(x) return np.arctanh(x) -def phi_second_derivative(x: float) -> float: +def phi_second_derivative(x: npt.NDArray[np.float64]) -> npt.NDArray[np.float64]: r"""Compute the second derivative of the DFT of the binary entropy. Parameters @@ -94,11 +104,10 @@ def phi_second_derivative(x: float) -> float: ------- An array corresponding to the image of the input array by the second derivative of the DFT of the binary entropy """ - x = np.asarray(x) return (1 - x**2) ** (-1) -def phi_inv(y: float, niter=3) -> float: +def phi_inv(y: npt.NDArray[np.float64]) -> npt.NDArray[np.float64]: r"""Compute the inverse of the DFT of the binary entropy via Halley's root finding algorithm Parameters @@ -114,34 +123,43 @@ def phi_inv(y: float, niter=3) -> float: An array corresponding to the image of the input array by the inverse of the DFT of the binary entropy. """ - y = np.asarray(y) + # We deal with the special case 0 and log(2) with masks + inverse = np.empty_like(y) + + inverse[y == np.log(2)] = 1 + inverse[y == 0] = 0 + + # We now compute the inverse for interior points + interior = (y > 0) & (y < np.log(2)) # Initial Guess of the inverse # Based on https://math.stackexchange.com/questions/3454390/find-the-approximation-of-the-inverse-of-binary-entropy-function - x = np.sqrt(1 - (1 - y / np.log(2)) ** (4 / 3)) + x = np.sqrt(1 - (1 - y[interior] / np.log(2)) ** (4 / 3)) # Iterative Halley's method for root finding with cubic convergence rate # See https://en.wikipedia.org/wiki/Halley's_method - for _ in range(niter): + for _ in range(3): - f = phi(x) - y + f = phi(x) - y[interior] df = phi_derivative(x) ddf = phi_second_derivative(x) x -= (f * df) / (df**2 - f * ddf / 2) - # We deal with the special case 0 and log(2) with where - return np.where(y > 0, np.where(y < np.log(2), x, 1), 0) + inverse[interior] = x + return inverse -def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: + +def mgl(mi_shares: npt.ArrayLike, group_order: int, base=2) -> npt.NDArray[np.float64]: r"""Upper bound the mutual information of a sensitive value in terms of the mutual information for each of its share. + Parameters ---------- - MI_shares : array_like, f64 + mi_shares : array_like, f64 Mutual information for each share. Array must be of shape ``(ns,nv)`` where - ``ns`` is the number of shares, ``nv`` the number of sensitive values. + ``n_shares`` is the number of shares, ``nv`` the number of sensitive values. group_order : int Order of the group in which the sensitive values are protected by masking. base : array_like, f64 @@ -152,37 +170,50 @@ def mrs_gerber_lemma(MI_shares: float, group_order: int, base=2) -> float: Upper bound on the mutual information for all nv sensitive values based on Mrs Gerber's Lemma. """ - MI_shares = np.asarray(MI_shares) + mi_shares = np.asarray(mi_shares, dtype=np.float64) + + if not (mi_shares >= 0).all(): + raise ValueError( + "Invalid inputs the mutual information on each share should be positive." + ) + if not (mi_shares <= np.log(group_order) / np.log(base)).all(): + raise ValueError( + "Invalid inputs the mutual information on each share should be less than the logarithm in base base of the group order." + ) # Check if the group order is a power of 2 (i.e. its bit expression contains a single 1) + # See Theorem 1 in the book chapter tosee the mgl reformulated using DFT as implemented here: + # Olivier Rioul, Julien Béguinot. The role of Mrs. Gerber’s Lemma for evaluating the information + # leakage of secret sharing schemes. Ioannis Kontoyiannis, Jason Klusowski, Cynthia Rush. Information + # Theory, Probability and Statistical Learning: A Festschrift in Honor of Andrew Barron, Springer, 2025. if (group_order & (group_order - 1)) == 0: - MI_share_bits = MI_shares * np.log2(base) - k = np.min(np.floor(MI_share_bits), axis=0) - clipped_MI_shares = np.clip(0, 1, MI_share_bits - k) * np.log(2) - MI = k * np.log(2) + phi(np.prod(phi_inv(clipped_MI_shares), axis=0)) + mi_share_bits = mi_shares * np.log2(base) + k = np.min(np.floor(mi_share_bits), axis=0) + clipped_mi_shares = np.clip(0, 1, mi_share_bits - k) * np.log(2) + mi = k * np.log(2) + phi(np.prod(phi_inv(clipped_mi_shares), axis=0)) # Otherwise, use a weaker MGL based on Pinsker/reverse Pinsker inequalities else: beta = group_order**2 * 4 ** (1 / group_order) # Detect which shares are below the noise amplification ratio. - below = 2 * MI_shares * np.log(base) < 1 + below = 2 * mi_shares * np.log(base) < 1 # Multipy only the terms bellow the noise amplification ratio - product = np.prod(np.where(below, 2 * MI_shares * np.log(base), 1), axis=0) / 4 + product = np.prod(np.where(below, 2 * mi_shares * np.log(base), 1), axis=0) / 4 # If there is no share bellow the noise amplification ratio we return the minimum instead - min_shares = np.min(MI_shares * np.log(base), axis=0) + min_shares = np.min(mi_shares * np.log(base), axis=0) # Depending on the case we return either the minimum or the product of shares below noise amplification P = np.where(~np.any(below, axis=0), min_shares, product) - MI_1 = np.log(1 + beta * P) - MI_2 = (1 / group_order + np.sqrt(P)) * np.log(1 + group_order * np.sqrt(P)) - min_MI = np.minimum(np.minimum(MI_1, MI_2), min_shares) + mi_1 = np.log(1 + beta * P) + mi_2 = (1 / group_order + np.sqrt(P)) * np.log1p(group_order * np.sqrt(P)) + min_mi = np.minimum(np.minimum(mi_1, mi_2), min_shares) # Depending on the case we return either the minimum or the amplification lemma - MI = np.where(~np.any(below, axis=0), min_shares, min_MI) + mi = np.where(~np.any(below, axis=0), min_shares, min_mi) # Conversion from nats to base 'base' - MI /= np.log(base) - return MI + mi /= np.log(base) + return mi diff --git a/tests/test_mgl.py b/tests/test_mgl.py index a4bef3bc..9b26e9f4 100644 --- a/tests/test_mgl.py +++ b/tests/test_mgl.py @@ -1,33 +1,54 @@ import pytest -from scalib.postprocessing import mrs_gerber_lemma +from scalib.postprocessing import mgl import numpy as np def test_mgl(): + rng = np.random.default_rng(seed=42) + err = 10**-15 - x = np.random.uniform(low=0, high=1) + x = rng.uniform(low=0, high=1) + + # Test the extreme case with group order of 2 + # Test with a single share + assert abs(mgl(x, 2, base=2) - x) <= err + + # Verify that a value error is raised here + with pytest.raises(ValueError): + mgl([x, 1.1], 2, base=2) + + with pytest.raises(ValueError): + mgl([1 + x, 1.1], 2, base=2) # Test with a single share - assert abs(mrs_gerber_lemma(x, 2**8, base=2) - x) <= err + assert abs(mgl(x, 2**8, base=2) - x) <= err # Test when a share leak more than a bit and another one less than a bit - assert abs(mrs_gerber_lemma([x, 1.1], 2**8, base=2) - x) <= err + assert abs(mgl([x, 1.1], 2**8, base=2) - x) <= err # Test when two shares leaks between 1 and 2 bits assert ( - abs( - mrs_gerber_lemma([1 + x, 1.1], 2**8, base=2) - - (1 + mrs_gerber_lemma([x, 0.1], 2**8, base=2)) - ) - <= err + abs(mgl([1 + x, 1.1], 2**8, base=2) - (1 + mgl([x, 0.1], 2**8, base=2))) <= err ) # Test when a single share is bellow noise amplification - assert abs(mrs_gerber_lemma(x, 2**8 - 1, base=2) - x) <= err + assert abs(mgl(x, 2**8 - 1, base=2) - x) <= err + + # Test when one share is bellow noise amplification threshold + assert abs(mgl([x, 5], 2**8 - 1, base=2) - x) <= err + + # Test when no share is bellow noise amplification threshold + assert abs(mgl([4, 5], 2**8 - 1, base=2) - 4) <= err + + # Test larger group order e.g. 3329 which the filed size in Kyber + q = 3329 + + # Test when a single share is bellow noise amplification + assert abs(mgl(x, q, base=2) - x) <= err # Test when one share is bellow noise amplification threshold - assert abs(mrs_gerber_lemma([x, 15], 2**8 - 1, base=2) - x) <= err + assert abs(mgl([x, 6], q, base=2) - x) <= err # Test when no share is bellow noise amplification threshold - assert abs(mrs_gerber_lemma([12, 15], 2**8 - 1, base=2) - 12) <= err + assert abs(mgl([5, 6], q, base=2) - 5) <= err From 1e0621e9046f4926757c258a9e41e07226bff039 Mon Sep 17 00:00:00 2001 From: beguinotj Date: Mon, 20 Oct 2025 14:09:17 +0200 Subject: [PATCH 12/12] slightly improved docs --- src/scalib/postprocessing/noise_amplification.py | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/scalib/postprocessing/noise_amplification.py b/src/scalib/postprocessing/noise_amplification.py index 14a02066..e70969c6 100644 --- a/src/scalib/postprocessing/noise_amplification.py +++ b/src/scalib/postprocessing/noise_amplification.py @@ -1,10 +1,14 @@ r"""Estimation of the mutual information between a sensitive value protected by masking and leakages in terms of the mutual information between each share and its corresponding leakages This function is usefull in the following setting. -You know that a sensitive value X is protected by masking so that it is shared into (S_0,...,S_d). -You can observe leakages Y_0,...,Y_d for each corresponding shares. -You have have estimated the leakages on each share via the mutual information I(S_i;Y_i). -Then, the mgl function provides an upper bound on the mutual information I(X; Y_0,...,Y_d). +You know that a sensitive value :math:`X` valued in group of order :math:`M` is protected by masking so that it is shared into :math:`(S_0,...,S_d)`. +You can observe leakages :math:`Y_0,...,Y_d` for each corresponding shares. +You have have estimated the leakages on each share via the mutual information :math:`I(S_i;Y_i)`. +Then, the mgl function provides an upper bound on the mutual information :math:`I(X; Y_0,...,Y_d)`: + +.. math:: + I(X; Y_0,\ldots,Y_d) \leq f_{\mathrm{MGL},M}( I(S_0;Y_0), \ldots, I(S_d;Y_d) ). + The obtained upper bound can then be used with other functions that provides security guarantees (such as success rate of an attack) in terms of mutual information. The `mgl` function takes as input the mutual information for each share separately (eventually for multiple sensitive values)