<jupyter_start><jupyter_text>*R*andom *S*ampling *plus* *F*ake *D*ata (RS+FD) ![image.png](attachment:image.png)  Some libraries<jupyter_code>import numpy as np

import pandas as pd
import matplotlib.pyplot as plt
import matplotlib
params = {'axes.titlesize':'15',
          'xtick.labelsize':'13',
          'ytick.labelsize':'13',
          'font.size':'15',
          'legend.fontsize':'medium',
          'lines.linewidth':'2',
          'font.weight':'normal',
          'lines.markersize':'10'
          }
matplotlib.rcParams.update(params)
plt.rcParams["mathtext.fontset"] = "cm"
%matplotlib inline

import sympy as sym
from mpl_toolkits.mplot3d import Axes3D<jupyter_output><empty_output><jupyter_text>Expectation of the estimator for RS+FD[GRR] - Eq. (4)<jupyter_code>n = sym.Symbol('n')
d = sym.Symbol('d')
k = sym.Symbol('k')
p = sym.Symbol('p')
q = sym.Symbol('q')
f = sym.Symbol('f')

E_n = (1/d) * (p*n*f + q*(n - n*f)) + (d-1)/(d*k) *n 
E_f = (d*E_n/(n*(p-q)) - (d-1+q*k)/(k*(p-q))) 
sym.simplify(E_f)<jupyter_output><empty_output><jupyter_text>Expectation of the estimator for RS+FD[OUE-z] - Eq. (8)<jupyter_code>n = sym.Symbol('n')
d = sym.Symbol('d')
k = sym.Symbol('k')
p = sym.Symbol('p')
q = sym.Symbol('q')
f = sym.Symbol('f')

E_n = (n/d) * (p*f + q*(1-f)) + (d-1)*n*q/d 
E_f = d/(n*(p-q)) *E_n - d*q/(p-q)
sym.simplify(E_f)<jupyter_output><empty_output><jupyter_text>Expectation of the estimator for RS+FD[OUE-r] -  Eq. (10)<jupyter_code>n = sym.Symbol('n')
d = sym.Symbol('d')
k = sym.Symbol('k')
p = sym.Symbol('p')
q = sym.Symbol('q')
f = sym.Symbol('f')

E_n = n/d * (f *(p-q) +q ) + (n*(d-1)/(d*k)) * (p-q+k*q)
E_f = d*E_n/(n*(p-q)) - ((p-q)*(d-1)+q*d*k)/(k*(p-q))
sym.simplify(E_f)<jupyter_output><empty_output><jupyter_text>RS+FD[GRR] Variance<jupyter_code>f = 0#sym.Symbol('f')

d = sym.Symbol('d')
k = sym.Symbol('k')
eps = sym.Symbol('\epsilon')
p_grr = sym.exp(eps)/(sym.exp(eps)+k-1)
q_grr = (1-p_grr)/(k-1)

sig_grr = (1/d) * (q_grr + f*(p_grr - q_grr) + (d-1)/k)
var_RSpFD_GRR= sym.simplify((d**2 * sig_grr * (1-sig_grr))/(p_grr-q_grr)**2)
var_RSpFD_GRR<jupyter_output><empty_output><jupyter_text>RS+FD[OUE-z] Variance<jupyter_code>f = 0#sym.Symbol('f')

d = sym.Symbol('d')
k = sym.Symbol('k')
eps = sym.Symbol('\epsilon')
p_oue = 0.5
q_oue = 1/(sym.exp(eps)+1)

sig_oue_z = (1/d) * (d*q_oue + f*(p_oue - q_oue))
var_RSpFD_OUE_z = sym.simplify((d**2 * sig_oue_z * (1-sig_oue_z))/(p_oue-q_oue)**2)
var_RSpFD_OUE_z<jupyter_output><empty_output><jupyter_text>RS+FD[OUE-r] Variance<jupyter_code>f = 0#sym.Symbol('f')

d = sym.Symbol('d')
k = sym.Symbol('k')
eps = sym.Symbol('\epsilon')
p_oue = 0.5
q_oue = 1/(sym.exp(eps)+1)

sig_oue_r = (1/d) * (f*(p_oue-q_oue) + q_oue + (d-1)*(p_oue/k + (k-1)*(q_oue/k)))
var_RSpFD_OUE_r = sym.simplify((d**2 * sig_oue_r * (1-sig_oue_r))/(p_oue-q_oue)**2)
var_RSpFD_OUE_r<jupyter_output><empty_output><jupyter_text>Analytical analysis: RS+FD[ADP]<jupyter_code>f = 0#sym.Symbol('f')

n = sym.Symbol('n')
d = sym.Symbol('d')
k = sym.Symbol('k')
eps = sym.Symbol('\epsilon')
p_grr = sym.exp(eps)/(sym.exp(eps)+k-1)
q_grr = (1-p_grr)/(k-1)

sig_grr = (1/d) * (q_grr + f*(p_grr - q_grr) + (d-1)/k)
var_RSpFD_GRR= sym.simplify((d**2 * sig_grr * (1-sig_grr))/(n*(p_grr-q_grr)**2))

p_oue = 0.5
q_oue = 1/(sym.exp(eps)+1)

sig_oue_z = (1/d) * (d*q_oue + f*(p_oue - q_oue))
var_RSpFD_OUE_z = sym.simplify((d**2 * sig_oue_z * (1-sig_oue_z))/(n*(p_oue-q_oue)**2))
resul = sym.simplify(var_RSpFD_GRR-var_RSpFD_OUE_z)
resul
def Variance_RSpFD_ADP(d, k, eps=np.log(3)):
    
    n= 10000 # number of users
    f=0
    
    #flipping probabilities with GRR
    p1 = np.exp(eps)/(np.exp(eps)+k-1)
    q1 = (1-p1)/(k-1)
    
    #flipping probabilities with OUE
    p2 = 1/2
    q2 = 1/(np.exp(eps) +1)
    
    delta1 = (q1 + f*(p1-q1) + (d-1)/k)
    d1md1 = delta1 * (1-delta1)
    
    delta2 = (d* q2 + f*(p2-q2))
    ineq = delta1*(d-delta1)*(p2-q2)**2 - delta2 * (d - delta2) * (p1 - q1) ** 2
    
    return -(2.0*d**2*k**2*(np.exp(eps) - 1)**2 + 2*(k + (d - 1)*(k + np.exp(eps) - 1))*(-d*k*(k + np.exp(eps) - 1) + k + (d - 1)*(k + np.exp(eps) - 1))*np.sinh(eps/2)**2)/(k**2*n*(np.exp(eps) - 1)**2*(np.cosh(eps) - 1))<jupyter_output><empty_output><jupyter_text>Generating Fig. 5 of the paper<jupyter_code>x = np.arange(2,11,1) #number of attributes
y = np.arange(2,21,1) #number of values per attribute
X, Y = np.meshgrid(x, y)
Z = Variance_RSpFD_ADP(X, Y)
fig = plt.figure(figsize = (10,5))
ax = plt.axes(projection='3d')

ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap='jet', edgecolor='lightgrey')
ax.set_xlabel('$d$', fontsize = 14)
ax.set_ylabel('$k_j$', fontsize = 14)
ax.set_yticks([5,10,15,20])
ax.set_zlabel('$VAR_1 - VAR_2$', fontsize = 14)
plt.subplots_adjust(top = 1, bottom = 0, right = 1, left = 0, hspace = 0, wspace = 0)
# plt.savefig('fig_surface_variance.png',dpi=500, bbox_inches = 'tight', pad_inches = 0.01)
plt.show()<jupyter_output><empty_output><jupyter_text>Analytical VS experimental<jupyter_code>def Var_RSpFD_GRR(n, d, k, p, q): #cf. Eq. (5)
    f_vi = 1/k  #uniform distribution
    
    sig_grr = (1/d) * (q + f_vi*(p-q) + (d-1)/k)
    
    var_grr = ((d**2 * sig_grr * (1-sig_grr))/(n*(p-q)**2))
    
    return var_grr

def Var_RSpFD_OUE_z(n, d, k, p, q): #cf. Eq. (9)
    f_vi = 1/k  #uniform distribution
    
    sig_oue_z = (1/d) * (d*q + f_vi * (p-q))
    
    var_oue_z = ((d**2 * sig_oue_z * (1-sig_oue_z))/(n*(p-q)**2))
    
    return var_oue_z

def Var_RSpFD_OUE_r(n, d, k, p, q): #cf. Eq. (11)
    f_vi = 1/k  #uniform distribution
    
    sig_oue_r = (1/d) * (f_vi * (p-q) + q + (d-1) * ((p/k) + ((k-1)/(k))*q) )
    
    var_oue_r = ((d**2 * sig_oue_r * (1-sig_oue_r))/(n*(p-q)**2))
    
    return var_oue_r<jupyter_output><empty_output><jupyter_text>Analytical Evaluation<jupyter_code># The values used for simulation here can be found in Subsection 4.1 - Setup of experiments

# Range of epsilon used in the experiments
lst_eps = [round(np.log(val),4) for val in range(2,8)] 

# Parameters of Synthetic Dataset 3. 
n = 50000 #number of users
d = 10 #number of attributes
beta = 1/d #sampling rate
K = [10 for _ in range(d)] #number of values per attribute is uniform, i.e., 10

# Initialize list of averaged variance values per epsilon
lst_var_RSpFD_GRR = []
lst_var_RSpFD_OUE_z = []
lst_var_RSpFD_OUE_r = []

for eps in lst_eps:
    eps_l = np.log((np.exp(eps) - 1 + beta)/beta) #amplified epsilon
    
    # Initialize list of variance values per attribute
    lst_var_att_grr = []
    lst_var_att_oue_z = []
    lst_var_att_our_r = []
    
    for k in K: # for each attribute
        
        p = np.exp(eps_l) / (np.exp(eps_l) + k - 1) # GRR
        q = (1 - p) / (k - 1) # GRR
        
        p_sue = 1/2 #OUE
        q_sue = np.exp(eps_l)/(np.exp(eps_l)+1) #OUE
        
        lst_var_att_grr.append(Var_RSpFD_GRR(n, d, k, p, q))
           
        lst_var_att_oue_z.append(Var_RSpFD_OUE_z(n, d, k, p_sue, q_sue))
        
        lst_var_att_our_r.append(Var_RSpFD_OUE_r(n, d, k, p_sue, q_sue))
           
    # We'll average the MSE over the number of attributes (cf. Eq. 13 of the paper)
    lst_var_RSpFD_GRR.append(np.mean(lst_var_att_grr))
    lst_var_RSpFD_OUE_z.append(np.mean(lst_var_att_oue_z))
    lst_var_RSpFD_OUE_r.append(np.mean(lst_var_att_our_r))<jupyter_output><empty_output><jupyter_text>Reading Experimental Evaluation results<jupyter_code># This requires executing the notebook 2_RS+FD_LDP_Synthetic3_50k_d10.ipynb

syn3_rs_fd_grr = pd.read_csv('results/synthetic3/final_results_rs_fd_grr.csv', header=None, index_col=0)
syn3_rs_fd_oue_z = pd.read_csv('results/synthetic3/final_results_rs_fd_oue_z.csv', header=None, index_col=0)
syn3_rs_fd_oue_r = pd.read_csv('results/synthetic3/final_results_rs_fd_oue_r.csv', header=None, index_col=0)
syn3_rs_fd_grr.shape, syn3_rs_fd_oue_z.shape, syn3_rs_fd_oue_r.shape #100 evaluations<jupyter_output><empty_output><jupyter_text>Plotting Curves<jupyter_code>fig, ax = plt.subplots(1,3, figsize=(20,4))
plt.subplots_adjust(wspace=0.28)

ax[0].plot(lst_var_RSpFD_GRR, label='RS+FD[GRR] - Ana')
ax[0].plot(syn3_rs_fd_grr.mean().values, label='RS+FD[GRR] - Exp')
ax[1].plot(lst_var_RSpFD_OUE_z, label='RS+FD[OUE-z] - Ana')
ax[1].plot(syn3_rs_fd_oue_z.mean().values, label='RS+FD[OUE-z] - Exp')
ax[2].plot(lst_var_RSpFD_OUE_r, label='RS+FD[OUE-r] - Ana')
ax[2].plot(syn3_rs_fd_oue_r.mean().values, label='RS+FD[OUE-r] - Exp')

for i in range(3):
    ax[i].set_xticks(range(len(lst_eps)))
    ax[i].set_xticklabels(lst_eps)
    ax[i].legend();
    ax[i].set_ylabel('Var')
    ax[i].set_xlabel('$\epsilon$')
    ax[i].set_yscale('log')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.19801064075873237, 'fraction_numerical': 0.03018736988202637, 'mean_word_length': 4.077510275983559, 'pattern_counts': {'":': 0, '<': 53, '<?xml version=': 0, '>': 53, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '640935', 'license': 'MIT', 'n_tokens_mistral': 4307, 'n_tokens_neox': 4121, 'n_words': 652, 'path': 'papers/[2]/1_RS+FD_Analysis.ipynb', 'repo_name': 'hharcolezi/ldp-protocols-mobility-cdrs', 'size': '8646'}