差分隐私实践入门

去标识

去标识与关联攻击

去标识（或匿名、假名）是指从数据集中删除标识信息的过程。例如有CSV数据集如下，第一行为表头：

1
2
3
4
Name,DOB,SSN,Zip,Workclass,Education,Education-Num,Marital Status,Occupation,Relationship,Race,Sex,Hours per week,Country,Target,Age,Capital Gain,Capital Loss
Karrie Trusslove,9/7/1967,732-14-6110,64152,State-gov,Bachelors,13,Never-married,Adm-clerical,Not-in-family,White,Male,40,United-States,<=50K,56,2174,0
Brandise Tripony,6/7/1988,150-19-2766,61523,Self-emp-not-inc,Bachelors,13,Married-civ-spouse,Exec-managerial,Husband,White,Male,13,United-States,<=50K,35,0,0
...

接下来导入并去标识：

1
2
3
4
5
6
7
import pandas
adult=pandas.read_csv("adult_with_pii.csv")
print(adult.head())

adult_data=adult.copy().drop(columns=['Name','SSN'])
adult_pii=adult[['Name','SSN','DOB','Zip']]
print(adult_data.head(1))

结果为：

1
2
3
4
5
6
7
8
9
10
11
12
               Name        DOB          SSN    Zip         Workclass  Education  ...  Hours per week        Country Target Age Capital Gain Capital Loss
0  Karrie Trusslove   9/7/1967  732-14-6110  64152         State-gov  Bachelors  ...              40  United-States  <=50K  56         2174            0        
1  Brandise Tripony   6/7/1988  150-19-2766  61523  Self-emp-not-inc  Bachelors  ...              13  United-States  <=50K  35            0            0        
2     Brenn McNeely   8/6/1991  725-59-9860  95668           Private    HS-grad  ...              40  United-States  <=50K  32            0            0        
3       Dorry Poter   4/6/2009  659-57-4974  25503           Private       11th  ...              40  United-States  <=50K  14            0            0        
4       Dick Honnan  9/16/1951  220-93-3811  75387           Private  Bachelors  ...              40           Cuba  <=50K  72            0            0        

[5 rows x 18 columns]
        DOB    Zip  Workclass  Education  Education-Num Marital Status  ... Hours per week        Country Target Age  Capital Gain Capital Loss
0  9/7/1967  64152  State-gov  Bachelors             13  Never-married  ...             40  United-States  <=50K  56          2174            0

[1 rows x 16 columns]

当尝试攻击的数据集与知道的一些辅助信息之间存在一些重叠列，可用这些列来实施一次简单的关联攻击：

1
2
karries_row=adult_pii[adult_pii['Name']=='Karrie Trusslove']
print(pandas.merge(karries_row,adult_data,left_on=['DOB','Zip'],right_on=['DOB','Zip']))

聚合

另一种防止隐私泄露的方法是只发布聚合数据，如只发布数据集平均年龄：

1
2
3
import pandas
adult=pandas.read_csv("adult_with_pii.csv")
print(adult['Age'].mean())

结果：

1
41.77250253355035

一般要将数据分组，给出给个分组的聚合统计信息：

1
print(adult[['Education','Age']].groupby('Education',as_index=False).mean().head(3))

结果：

1
2
3
4
  Education        Age
0      10th  42.032154
1      11th  42.057021
2      12th  41.879908

如果某分组个体较少，聚合统计结果将难以达到隐私保护目的。此时采用差分攻击，对数据集中某个大分组执行两次求和问询。第一次对整个数据集进行问询，第二次除一条记录外所有记录进行问询。

1
print(adult['Age'].sum()-adult[adult['Name']!='Karrie Trusslove']['Age'].sum())

结果：

1
56

$k$-匿名性

对于特定的$k$，若对于任意记录$r_1\in D$，存在至少$k-1$条其他的记录$r_2,\cdots,r_k\in D$，使得$\displaystyle\prod_{\operatorname{qi}(D)}r_1=\prod_{\operatorname{qi}(D)}r_2,\cdots,\prod_{\operatorname{qi}(D)}r_k$，其中$\operatorname{qi}(D)$是$D$的准标识，$\displaystyle\prod_{\operatorname{qi}(D)}r$标识包含准标识的列$r$，则称数据集$D$满足$k$-匿名性。

换句话说，一部分辅助数据不应该“过多地”缩小个题所属记录的可能范围，目的是保证每个个题都能“融入人群”。

把数据集按照数据集各列中的特定子集分组，即按照准标识分组，是每个分组中的个体都拥有相同的准标识。若数据集中每个个体所属分组大小都至少为$k$，则称此数据集满足$k$-匿名性。虽攻击者仍可将攻击范围缩小至特定分组中，但无法进一步确定分组中哪个个体才是攻击目标。

此时求$k$：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
import pandas,numpy
import matplotlib.pyplot
matplotlib.pyplot.style.use('seaborn-whitegrid')
raw_data={
    'first_name':['Jason','Molly','Tina','Jake','Amy'], 
    'last_name':['Miller','Jacobson','Ali','Milner','Cooze'], 
    'age':[42,52,36,24,73], 
    'preTestScore':[4,24,31,2,3],
    'postTestScore':[25,94,57,62,70]
}
#df = pd.DataFrame(raw_data, columns = ['first_name', 'last_name', 'age', 'preTestScore', 'postTestScore'])
df=pandas.DataFrame(raw_data,columns=['age','preTestScore','postTestScore'])
print(df)
def isKAnonymized(df,k):
    for index,row in df.iterrows():
        query=' & '.join([f'{col}=={row[col]}' for col in df.columns])
        rows=df.query(query)
        if rows.shape[0]<k:
            return False
    return True
print(isKAnonymized(df,1))
print(isKAnonymized(df,2))

结果如下，说明满足$1$-匿名性：

1
2
3
4
5
6
7
8
   age  preTestScore  postTestScore
0   42             4             25
1   52            24             94
2   36            31             57
3   24             2             62
4   73             3             70
True
False

可通过泛化数据的方式对数据集进行修改，使其满足特定取值$k$下的$k$-匿名性。泛化指的是将数据修改为不那么特殊的数据，使其更可能与数据集中其他个体的数据相匹配。例如精确到个位的年龄可通过四舍五入泛化为精确到十位；将邮政编码最右侧数字替换为$0$来泛化邮政编码等。

1
2
3
4
5
6
7
8
9
def generalize(df,depths): #depths称为查找表 存储每一列要用0替换多少位数字
    return df.apply(lambda x:x.apply(lambda y:int(int(y/(10**depths[x.name]))*(10**depths[x.name]))))
depths={
    'age':1,
    'preTestScore':1,
    'postTestScore':1
}
df2=generalize(df,depths) #一层泛化 四舍五入到十位
print(df2)

结果为：

1
2
3
4
5
6
   age  preTestScore  postTestScore
0   40             0             20
1   50            20             90
2   30            30             50
3   20             0             60
4   70             0             70

尽管还是无法满足$2$-匿名性，但二层泛化会删除所有数据。在拥有更多个体的数据集中，通常需要更少的泛化处理即可使数据集满足所需$k$取值下的$k$-匿名性。

1
2
3
4
5
6
7
8
9
10
11
12
adult_data=pandas.read_csv("adult_with_pii.csv")
print(adult_data.head())
df=adult_data[['Age','Education-Num']]
df.columns=['age','edu']
print(isKAnonymized(df.head(100),1))
print(isKAnonymized(df.head(100),2))
depths={
    'age':1,
    'edu':1
}
df2=generalize(df.head(100),depths)
print(isKAnonymized(df2,2))

结果：

1
2
3
4
5
6
7
8
9
10
11
               Name        DOB          SSN    Zip         Workclass  Education  ...  Hours per week        Country Target Age Capital Gain Capital Loss
0  Karrie Trusslove   9/7/1967  732-14-6110  64152         State-gov  Bachelors  ...              40  United-States  <=50K  56         2174            0        
1  Brandise Tripony   6/7/1988  150-19-2766  61523  Self-emp-not-inc  Bachelors  ...              13  United-States  <=50K  35            0            0        
2     Brenn McNeely   8/6/1991  725-59-9860  95668           Private    HS-grad  ...              40  United-States  <=50K  32            0            0        
3       Dorry Poter   4/6/2009  659-57-4974  25503           Private       11th  ...              40  United-States  <=50K  14            0            0        
4       Dick Honnan  9/16/1951  220-93-3811  75387           Private  Bachelors  ...              40           Cuba  <=50K  72            0            0        

[5 rows x 18 columns]
True
False
False

可见还是不满足$2$-匿名性，根本原因是数据集包含异常值，包含一些与其他个体差异非常大的个体，即使泛化处理后也很难使这些异常个体融入任何分组中。在这里，非常低和非常高的年龄在数据集中过于稀少，进一步泛化每一行数据又会过分泛化年龄在$[20,40]$范围内的数据，损害结果数据可用性。实际上，实现满足$k$-匿名性的最优泛化方法被证明是个NP-困难问题。

异常值问题的一个简单解决方案是将数据集中每个个体年龄限制在一个特定范围内，从而完全消除数据集中的异常值。但这种方法会损害数据可用性，因为要用假年龄值替代真实年龄值。在这里将所有年龄取值都限制在$60$岁及以下，且不再考虑受教育年数这一列数据，都替换成一个非常大的值。

1
2
3
4
5
6
7
8
depths={
    'age':1,
    'edu':1
}
dfp=df.clip(upper=numpy.array([60,100000000]),axis='columns')
dfp=dfp.clip(lower=numpy.array([10,3]),axis='columns')
df2=generalize(dfp.head(500),depths)
print(isKAnonymized(df2,7))

结果可见满足$7$-匿名性：

1
True

差分隐私

满足差分隐私的函数称为机制，若对所有邻近数据集$x$和$x’$，以及所有可能的输出$S$，机制$F$均满足$\dfrac{\operatorname{Pr}(F(x)=S)}{\operatorname{Pr}(F(x’)=S)}\leqslant e^{\epsilon}$，则称机制$F$满足差分隐私。

差分隐私不是数据具有的属性，而是算法具有的属性，我们可证明一个算法满足差分隐私。若想证明一个数据集满足差分隐私，需要证明产生此数据集的算法满足差分隐私。

若两个数据集中只有一个个体的数据项不同，则认为这两个数据集是临近数据集。$\operatorname{F}$是个随机函数，给定相同输入一般也包含多个可能的输出，概率分布一般不是点分布。$F$引入的随机性应该足够大，使得观察$F$的输出无法判断输入是$x$还是$x’$。假设某数据在$x$中但不在$x’$中，如果攻击者无法确定输入是$x$还是$x’$，则无法判断输入是否包含该数据，更不用说判断该数据是什么了。

参数$\epsilon$称为隐私参数（或隐私预算），$\epsilon$较小时意味着$F$需要为相似的输入提供非常相似的输出，因此提更高等级的隐私性。实际中一般设为约等于$1$或更小的值，大于$10$的值一般意味着大概率无法提供足够的隐私性。

差分隐私一般用于回复特定问询，如数据集中有多少个体年龄大于等于$40$岁：

1
2
3
import pandas,numpy
adult=pandas.read_csv("adult_with_pii.csv")
print(adult[adult['Age']>=40].shape[0])

结果：

1
17450

使这个问询满足差分隐私的最简单方法实在回复结果上增加随机噪声，既需要增加足够大的噪声，又不能加得太多，否则问询结果就无意义了。此时一些基础机制具体描述了应增加何种类型的噪声，以及噪声量应该有多大。

根据拉普拉斯机制，对于可以输出一个数值型结果的函数$f(x)$，定义$F(x)=f(x)+\operatorname{Lap}\left(\dfrac{s}{\epsilon}\right)$，其中$s$是$f$的敏感度，$\operatorname{Lap}(S)$表示以均值为$0$、放缩系数为$S$的拉普拉斯分布采样。

函数$f$的敏感度是指当输入由数据集$x$变化为邻近数据集$x’$后，$f$的输出变化量。计数问询的敏感度总为$1$，因为当问询数据集中满足特定属性的数据量时，如果只修改数据集中的一个数据项，则问询的输出变化量最多为$1$。

1
2
3
sensitivity=1
epsilon=0.1
print(adult[adult['Age']>=40].shape[0]+numpy.random.laplace(loc=0,scale=sensitivity/epsilon))

结果不唯一：

1
17432.36769513577

此时构造一个恶意的计数问询，确定某个人的收入是否大于$50000$$。

1
2
3
4
5
6
7
karries_row=adult[adult['Name']=='Karrie Trusslove']
print(karries_row[karries_row['Target']=='<=50K'].shape[0])

sensitivity=1
epsilon=0.1
karries_row=adult[adult['Name']=='Karrie Trusslove']
print(karries_row[karries_row['Target']=='<=50K'].shape[0]+numpy.random.laplace(loc=0,scale=sensitivity/epsilon))

结果不唯一：

1
2
1
13.558463386047062

差分隐私的性质

差分隐私的一个性质是串行组合性。如果$F_1(x)$满足$\epsilon_1$-差分隐私，$F_2(x)$满足$\epsilon_2$-差分隐私，则同时发布两个结果的机制$G(x)=(F_1(x),F_2(x))$满足$(\epsilon_1+\epsilon_2)$-差分隐私。

1
2
3
4
5
6
7
8
9
10
11
12
import numpy
epsilon1=1
epsilon2=1
epsilon_total=2
def F1():
    return numpy.random.laplace(loc=0,scale=1/epsilon1)
def F2():
    return numpy.random.laplace(loc=0,scale=1/epsilon2)
def F3():
    return numpy.random.laplace(loc=0,scale=1/epsilon_total)
def F_combined(): #满足2-差分隐私
    return (F1()+F2())/2

串行组合性给出了多次发布后总$\epsilon$的上界，实际发布的积累隐私消耗量可能会更低一些，所以F_combined噪声值分布图可能比F3更平一些。

差分隐私的另一个性质为并行组合性。如果$F(x)$满足$\epsilon$-差分隐私，将数据集$X$切分为$k$个互不相交的子数据块$\displaystyle\bigcup_{i=1}^kx_i=X$，则发布所有结果$F(x_1),F(x_2),\cdots,F(x_k)$的机制满足$\epsilon$-差分隐私。

直方图应用并行组合性：

1
2
3
4
5
6
7
8
import pandas,numpy
adult=pandas.read_csv("adult_with_pii.csv")
print(adult['Education'].value_counts().to_frame().head(5))

epsilon = 1
f=lambda x:x+numpy.random.laplace(loc=0,scale=1/epsilon)
s=adult['Education'].value_counts().apply(f)
print(s.to_frame().head(5))

结果：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
              count
Education
HS-grad       10501
Some-college   7291
Bachelors      5355
Masters        1723
Assoc-voc      1383
                     count
Education
HS-grad       10501.113144
Some-college   7290.343125
Bachelors      5355.624063
Masters        1725.033860
Assoc-voc      1382.385230

列联表（或交叉列表、交叉表）应用并行组合性：

1
2
3
4
5
6
7
8
import pandas,numpy
adult=pandas.read_csv("adult_with_pii.csv")
print(pandas.crosstab(adult['Education'],adult['Sex']).head(5))

epsilon=1
ct=pandas.crosstab(adult['Education'],adult['Sex'])
f=lambda x:x+numpy.random.laplace(loc=0,scale=1/epsilon)
print(ct.applymap(f).head(5))

结果：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
Sex        Female  Male
Education
10th          295   638
11th          432   743
12th          144   289
1st-4th        46   122
5th-6th        84   249
Sex            Female        Male
Education
10th       295.088175  639.521246
11th       434.281535  742.165123
12th       144.187342  288.693666
1st-4th     46.075385  121.528562
5th-6th     84.914934  247.715513

但若为列联表每增加一个属性后，将数据集分割称更多数据块，每个数据块包含的数据行数会变少，每个计数结果都会变小。真实计数值的缩小对差分隐私计数结果的准确性造成巨大影响。较大的计数结果表示一个强信号，不太可能被一个相对较弱的噪声所淹没；反之较小的计数结果表示一个弱信号，强度甚至可能接近噪声值，就无法从结果中获得任何有用信息了。

另一个差分隐私的性质为后处理性，即不可能通过某种方式对差分隐私保护下的数据进行后处理，来降低差分隐私的隐私保护程度。若$F(x)$满足$\epsilon$-差分隐私，则对于任意确定或随机函数$g$，$g(F(x))$也满足$\epsilon$-差分隐私。还可以对计算结果进行后处理，以降低噪声量、改善输出成果。

后处理性的另一个含义是，差分隐私可抵抗基于辅助数据的隐私攻击方式。如$g$可能包含关于数据集元素的辅助数据，$g$期望利用该信息实施关联攻击，但后处理性意味着攻击效果会被隐私参数$\epsilon$约束。

敏感度

在拉普拉斯机制中，使问询满足差分隐私所需噪声量取决于问询的敏感度。给定一个将数据集$\mathcal D$映射为实数的函数$f:\mathcal D\rightarrow\mathbb R$，$f$的全局敏感度为$\displaystyle\operatorname{GS}(f)=\max_{x,x’:d(x,x’)\leqslant1}\left|f(x)-f(x’)\right|$。即对于任意两个邻近数据集$x$和$x’$，$f(x)$和$f(x’)$最多相差$\operatorname{GS}(f)$。其中$d(x,x’)$表示两个数据集$x$和$x’$之间的距离。若两个数据集之间距离小于等于$1$，则称这两个数据集是邻近集。

对于距离$d(x,x’)$，这里使用两个数据集的对称差$d(x,x’)=|x-x’\cup x’-x|$。即如果$x’$是通过在$x$中添加或删除一行来构造的，则$d(x,x’)=1$；如果$x’$是通过在$x$中修改一行来构造的，则$d(x,x’)=2$。

当$x$变化$1$，$f(x)$变化的值即为全局敏感度。当$f(x)$的变化取决于$x$，则称全局敏感度是无界的。

计数问询计算数据集中满足特定属性的行数，敏感度总等于$1$。

1
2
3
4
5
6
import pandas
adult=pandas.read_csv("adult_with_pii.csv")
print(adult.shape[0])
print(adult[adult['Education-Num']>10].shape[0])
print(adult[adult['Education-Num']<=10].shape[0])
print(adult[adult['Name']=='Joe Near'].shape[0])

结果：

1
2
3
4
32563
10517
22046
0

求和问询计算数据集行中的属性值总和。想数据集添加一个新数据会使问询增加一个新个体的年龄，因此敏感度取决于添加的具体数据是什么。此时敏感度很难找到一个确切的数，例如敏感度设为$125$，有可能年龄会有$126$岁的人。对于收入等其他属性来说，更难找到一个合理的上界。当待求和的属性值不存在上下界时，称求和问询具有无界敏感度。存在上下界时求和问询敏感度等于上下界的差。

1
print(adult[adult['Education-Num']>10]['Age'].sum())

结果：

1
441431

均值问询计算指定列属性值的平均值，可拆分成求和问询除以计数问询，两者敏感度都可计算。可以分别计算使用拉普拉斯机制的两个问询的噪声回复并做触发，即可得到差分隐私均值，并应用串行组合性计算两个问询的总隐私消耗量。

1
2
print(adult[adult['Education-Num']>10]['Age'].mean())
print(adult[adult['Education-Num']>10]['Age'].sum()/adult[adult['Education-Num']>10]['Age'].shape[0])

结果：

1
2
41.973091185699346
41.973091185699346

差分隐私的拉普拉斯机制无法直接应用于无界敏感度问询，通常用裁剪技术将此类问询转换为等价的有界敏感度问询，即强制设置属性值的上下界。例如$125$岁以上的年龄被裁剪到恰好为$125$岁。对裁剪数据执行求和问询的敏感度等于上下界的差。

1
print(adult['Age'].clip(lower=0,upper=125).sum()) #敏感度125

结果：

1
1360238

裁剪边界的上下界越接近，敏感度越低，差分隐私所需噪声就越小。过分裁剪会从数据中移除很多信息，给准确度带来损失，影响超过降低敏感度改善噪声所带来的准确度提升效果。

一般来说，要将裁剪边界设置为尽可能100%保留数据集所有信息，然而某些领域很难做到这一点。通过直方图寻找数据集的上界不满足差分隐私，因为边界可能会泄露数据的一些相关信息。一般根据数据集先天满足的一些属性来确定裁剪边界，即无需查看数据本身即可获得此属性。

也可先将敏感度下界设为$0$，逐渐增加上界，直至问询输出不再变化。问询结果应当开始逐渐增大，之后趋于平缓，最后噪声淹没信号，波动非常剧烈。

近似差分隐私

近似差分隐私（或$(\epsilon,\delta)$-差分隐私）定义为$\operatorname{Pr}(F(x)=S)\leqslant e^{\epsilon}\operatorname{Pr}(F(x’)=S)+\delta$。隐私参数$\delta$表示不满足此近似差分隐私定义的失败概率。此时有$1-\delta$的概率获得等价于纯粹差分隐私的隐私保护程度，有$\delta$概率不满足隐私参数为$\epsilon$的纯粹差分隐私。要求$\delta$很小，通常$\delta\leqslant\dfrac{1}{(\#\mathcal D)^2}$。

近似差分隐私同样满足串并行组合性和后处理性。但对于串行组合性，若$F_1(x)$和$F_2(x)$分别满足$(\epsilon_1,\delta_1)$-差分隐私和$(\epsilon_2,\delta_2)$-差分隐私，则发布两个结果的机制$G(x)=(F_1(x),F_2(x))$满足$(\epsilon_1+\epsilon_2,\delta_1+\delta_2)$-差分隐私。

高斯机制无法满足纯粹$\epsilon$-差分隐私，但可以满足$(\epsilon,\delta)$-差分隐私。对于返回数值的函数$f(x)$，应用$F(x)=f(x)+\mathcal N\left(0,\sigma^2\right)$，得到满足$(\epsilon,\delta)$-差分隐私的$F(x)$。这里$s$是$f$的敏感度，$\mathcal N\left(0,\sigma^2\right)$表示均值为$0$，方差为$\sigma^2$的高斯分布采样结果，其中$\sigma^2=\dfrac{2s^2\ln\left(\dfrac{1.25}{\delta}\right)}{\epsilon^2}$。

给定长度为$k$的向量$\overrightarrow V$，其L1范数定义为$\displaystyle\left|\left|\overrightarrow V\right|\right|_1=\sum_{i=1}^k\overrightarrow{V_i}$，二维空间中$\left|\left|\overrightarrow{V_1}-\overrightarrow{V_2}\right|\right|_1$就是他们的曼哈顿距离。L2范数定义为$\left|\left|\overrightarrow V\right|\right|_2=\sqrt{\displaystyle\sum_{i=1}^k\overrightarrow{V_i^2}}$，二维空间中$\left|\left|\overrightarrow{V_1}-\overrightarrow{V_2}\right|\right|_2$就是他们的欧几里得距离。有$\left|\left|\overrightarrow{V_1}-\overrightarrow{V_2}\right|\right|_1\geqslant\left|\left|\overrightarrow{V_1}-\overrightarrow{V_2}\right|\right|_2$。

拉普拉斯机制和高斯机制都可以扩展到形式为$f:\mathcal D\rightarrow{\mathbb R}^k$的向量值函数。向量值拉普拉斯机制需要使用L1敏感度，而向量值高斯机制可以使用L1或L2敏感度。

向量值函数$f$的L1敏感度为$\displaystyle\operatorname{GS}_1=\max_{d(x,x’)\leqslant1}\left|\left|f(x)-f\left(x’\right)\right|\right|_1$，若$f$返回一个长度为$k$的向量，且向量中各个元素的敏感度均为$1$，则$f$的L1敏感度为$k$。向量值函数$f$的L2敏感度为$\displaystyle\operatorname{GS}_2=\max_{d(x,x’)\leqslant1}\left|\left|f(x)-f\left(x’\right)\right|\right|_2$，若$f$返回一个长度为$k$的向量，且向量中各个元素的敏感度均为$1$，则$f$的L2敏感度为$\sqrt k$。

向量拉普拉斯机制和向量高斯机制发布的都是$f(x)+\left(Y_1,Y_2,\cdots,Y_k\right)$，$Y_i$是独立同分布噪声。

高级组合定理通常用$k$-折叠适应性组合机制来描述，后者指的是将一系列机制$m_1,m_2,\cdots,m_k$组合起来，这些机制满足：

  • 适应性：可根据所有前述机制$m_1,m_2,\cdots,m_{i-1}$的输出来选择下一个机制$m_i$。
  • 组合性：每个机制$m_i$的输入既包括隐私数据集，也包括前述机制的所有输出。

如果$k$-折叠适应性组合$m_1,m_2,\cdots,m_k$中的每个机制$m_i$都满足$(\epsilon,\delta)$-差分隐私，则$\forall\delta’\geqslant0$，整个$k$-折叠适应性组合满足$\left(\epsilon’,k\delta+\delta’\right)$-差分隐私，其中$\epsilon’=2\epsilon\sqrt{2k\ln\left(\dfrac{1}{\delta’}\right)}$。当待组合的机制满足纯粹$\epsilon$-差分隐私时，有$\delta=0$。

对于相同机制，应用高级组合性得到的$\epsilon’$下届远低于应用串行组合性得到的下界，意味着串行组合性得到的隐私消耗量下界是宽松的，高级组合性得到的下界稍显紧致一些。事实证明当$k<70$时，标准的串行组合性比高级组合性得到的总隐私消耗量更小。当$k>100$时，高级组合性可显露出巨大的优势。

局部敏感度

全局敏感度定义考察的是任意两个邻近数据集，但由于将在实际数据集上执行差分隐私机制，应该只需要考虑此数据集的邻近数据集。由此得到函数$f:\mathcal D\rightarrow\mathbb R$在$x:\mathcal D$的局部敏感度定义为$\displaystyle\operatorname{LS}(f,x)=\max_{x’:d(x,x’)\leqslant1}\left|f(x)-f\left(x’\right)\right|$。

建议—测试—发布框架先询问数据分析者函数的建议局部敏感度上界，然后执行满足差分隐私的测试，检验所问询的数据集是否远离了局部敏感度高于建议边界的数据集。若测试通过，该框架发布噪声结果，并将噪声量校准到建议的边界。建议—测试—发布框架具体如下，满足$(\epsilon,\delta)$-差分隐私：

  • 建议一个局部敏感度的目标边界$b$。
  • 若$D(f,x,b)+\operatorname{Lap}\left(\dfrac{1}{\epsilon}\right)<\dfrac{\ln\left(\dfrac{2}{\delta}\right)}{2\epsilon}$，返回$\bot$。
  • 否则返回$f(x)+\operatorname{Lap}\left(\dfrac{b}{\epsilon}\right)$。

这里用$A(f,x,k)$表示通过从数据集$x$执行$k$步可得到$f$的最大局部敏感度，即$\displaystyle A(f,x,k)=\max_{y:d(x,y)\leqslant k}\operatorname{LS}(f,y)$。然后定义$k$距离局部敏感度表示需要多少步才能实现比给定上界$b$更大的局部敏感度，即$D(f,x,b)=\underset{k}{\operatorname{argmin}}(A(f,x,k)>b)$，$D(f,x,b)$的全局敏感度为$1$。

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
import pandas,numpy
def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def ls_at_distance(df,u,k):
    return numpy.abs(u/(len(df)-k+1))
def dist_to_high_ls(df,u,b):
    k=0
    while ls_at_distance(df,u,k)<b:
        k+=1
    return k
def ptr_avg(df,u,b,epsilon,delta,logging=False):
    df_clipped=df.clip(upper=u)
    k=dist_to_high_ls(df_clipped,u,b)
    noisy_distance=laplace_mech(k,1,epsilon)
    threshold=numpy.log(2/delta)/(2*epsilon)
    if logging:
        print(f"噪声距离：{noisy_distance}，门限值：{threshold}")
    if noisy_distance>=threshold:
        return laplace_mech(df_clipped.mean(),b,epsilon)
    else:
        return None
adult=pandas.read_csv("adult_with_pii.csv")
df=adult['Age']
u=100 #年龄上界
epsilon=1
delta=1/(len(df)**2)
b=0.005 #建议敏感度
print(ptr_avg(df,u,b,epsilon,delta,logging=True))

结果为：

1
2
噪声距离：12563.298110920712，门限值：10.737505543744179
41.79231612235209

采样—聚合框架：对任何函数$f:\mathcal D\rightarrow\mathbb R$，令裁剪上下界分别为$u$和$l$，则下述框架满足$\epsilon$-框架：

  • 将数据集$X\in\mathcal D$拆分成不相交的数据块$x_1,x_2,\cdots,x_k$。
  • 计算每个数据块的裁剪回复值$a_i=\max(l,\min(u,f(x_i)))$。
  • 计算平均回复值并增加噪声$\displaystyle A=\left(\dfrac{1}{k}\sum_{i=1}^ka_i\right)+\operatorname{Lap}\left(\dfrac{u-l}{k\epsilon}\right)$。

分块数量$k$越大，噪声均值敏感度越小，噪声量越小，但回复值$f(x_i)$偏离正确回复值$f(X)$越远。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import pandas,numpy
def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def f(df):
    return df.mean()
def saa_avg_age(k,epsilon,logging=False):
    df=adult['Age']
    chunk_size=int(numpy.ceil(df.shape[0]/k)) #计算每个数据块应包含的行数
    if logging:
        print(f'Chunk size: {chunk_size}')
    xs=[df[i:i+chunk_size] for i in range(0,df.shape[0],chunk_size)] #将df拆分为数据块
    answers=[f(x_i) for x_i in xs] #在每个x_i上执行f 并裁剪输出值
    u=80
    l=20
    clipped_answers=numpy.clip(answers,l,u)
    noisy_mean=laplace_mech(numpy.mean(clipped_answers),(u-l)/k,epsilon) #计算输出均值 并增加噪声
    return noisy_mean
adult=pandas.read_csv("adult_with_pii.csv")
print(saa_avg_age(600,1,logging=True))

结果：

1
2
Chunk size: 55
41.74177911583807

差分隐私变体

散度是度量两种概率分布差异的方法，最大散度是Kullback-Leibler散度的最坏情况，两个概率分布$Y$和$Z$的最大散度定义为$\displaystyle D_{\infty}(Y\mid\mid Z)=\max_{S\subseteq\operatorname{Supp}(Y)}\ln\dfrac{\operatorname{Pr}(Y\in S)}{\operatorname{Pr}(Z\in S)}$。如果证明$D_{\infty}(F(x)\mid\mid F(x’))\leqslant\epsilon$则$F$满足$\epsilon$-差分隐私。

概率分布$P$和$Q$的$\alpha$阶Rényi散度定义为$D_{\infty}(P\mid\mid Q)=\dfrac{1}{\alpha-1}\ln \left(E_{x\sim Q}\left(\dfrac{P(x)}{Q(x)}\right)^{\alpha}\right)$，其中$P(x)$和$Q(x)$为在点$x$处的概率密度。

Rényi差分隐私RDP：若对于所有邻近数据集$x$和$x’$，随机机制$F$满足$D_{\alpha}(F(x)\mid\mid F(x’))\leqslant\overline{\epsilon}$，则称此机制$F$满足$(\alpha,\epsilon)$-RDP，这里$\overline{\epsilon}$和之前的$\epsilon$没啥关系。若$F$满足$(\alpha,\overline{\epsilon})$-RDP，则对任意$\delta>0$，$F$满足$(\epsilon,\delta)$-差分隐私，其中$\epsilon=\overline{\epsilon}+\dfrac{\ln\left(\dfrac{1}{\delta}\right)}{\alpha-1}$。实现Rényi差分隐私的基本机制的高斯机制。

对于一个L2敏感度为$\Delta f$的函数$f:\mathcal D\rightarrow\mathbb{R}^k$，可构造$(\alpha,\overline{\epsilon})$-RDP机制：$F(x)=f(x)+\mathcal N(\sigma^2),\sigma^2=\dfrac{\Delta f^2\alpha}{2\overline{\epsilon}}$。

1
2
3
def gaussian_mech_RDP_vec(vec,sensitivity,alpha,epsilon_bar):
    sigma=numpy.sqrt((sensitivity**2*alpha)/(2*epsilon_bar))
    return [v+numpy.random.normal(loc=0,scale=sigma) for v in vec]

Rényi差分隐私的串行组合性描述为：如果$F_1$和$F_2$分别满足$(\alpha,\overline{\epsilon_1})$-RDP和$(\alpha,\overline{\epsilon_2})$-RDP，则组合机制满足$(\alpha,\overline{\epsilon_1}+\overline{\epsilon_2})$-RDP。当给定噪声等级$\sigma^2$时，用Rényi差分隐私的串行组合性来限制重复应用高斯机制的隐私消耗量，再将隐私定义转换为$(\epsilon,\delta)$-差分隐私，总隐私消耗量比直接在$(\epsilon,\delta)$中应用串行组合定理或高级组合定理得到的总隐私消耗量要低得多。

零集中差分隐私zCDP：如果对于所有邻近数据集$x$和$x’$，以及所有的$\alpha\in(1,\infty)$，随机机制$F$满足$D_{\alpha}(F(x)\mid\mid F(x’))\leqslant\rho\alpha$，则称$F$满足$\rho$-zCDP。若$F$满足$\rho$-zCDP，则对于任意给定的$\delta>0$，$F$满足$(\epsilon,\delta)$-差分隐私，其中$\epsilon=\rho+2\sqrt{\rho\ln\left(\dfrac{1}{\delta}\right)}$。零集中差分隐私也是根据Rényi散度定义的差分隐私变体。

对于一个L2敏感度为$\Delta f$的函数$f:\mathcal D\rightarrow\mathbb{R}^k$，可构造$\rho$-zCDP机制：$F(x)=f(x)+\mathcal N(\sigma^2),\sigma^2=\dfrac{\Delta f^2}{2\rho}$。

1
2
3
def gaussian_mech_zCDP_vec(vec,sensitivity,rho):
    sigma=numpy.sqrt((sensitivity**2)/(2*rho))
    return [v+numpy.random.normal(loc=0,scale=sigma) for v in vec]

零集中差分隐私的串行组合性描述为：若$F_1$和$F_2$分别满足$\rho_1$-zCDP和$\rho_2$-zCDP，则组合机制满足$(\rho_1+\rho_2)$-zCDP。

Rényi差分隐私和零集中差分隐私为同一算法提供了更加进制的隐私消耗量上界，因此使用高斯机制构造迭代算法时有限考虑使用这些差分隐私变体。

指数机制

对于非数值型回复，无法直接在结果上增加噪声，还要保证回复过程满足差分隐私。指数机制可以从备选回复集合中选出最佳回复的同时，保证回复过程满足差分隐私。分析者需要定义一个备选回复集合$\mathcal R$，并指定一个全局敏感度为$\Delta u$评分函数$u:\mathcal D\times\mathcal R\rightarrow\mathbb R$，输出备选回复集合中每个回复的分数，分数最高的就是最佳回复$r\in\mathcal R$，各个回复的输出概率与$\exp\left(\dfrac{\epsilon u(x,r)}{2\Delta u}\right)$成正比。指数机制通过返回分数近似最大的回复来实现差分隐私保护，即返回的结果所对应的分数可能不是备选回复集合中分数最高的那个结果。

无论$\mathcal R$中包含多少个备选输出，指数机制的总隐私消耗量仍然为$\epsilon$。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
import pandas,numpy
adult=pandas.read_csv("adult_with_pii.csv")
def score(data,option):
    return data.value_counts()[option]/1000
print(score(adult['Marital Status'], 'Never-married'))

def exponential(x,R,u,sensitivity,epsilon):
    scores=[u(x,r) for r in R] #计算R中每个回复的分数
    probabilities=[numpy.exp(epsilon*score/(2*sensitivity)) for score in scores] #根据分数计算每个回复的输出概率
    probabilities=probabilities/numpy.linalg.norm(probabilities,ord=1) #对概率进行归一化处理 使概率和等于1
    return numpy.random.choice(R,1,p=probabilities)[0] #根据概率分布选择回复结果
options=adult['Marital Status'].unique()
print(exponential(adult['Marital Status'],options,score,1,1))

r=[exponential(adult['Marital Status'],options,score,1,1) for i in range(200)]
print(pandas.Series(r).value_counts())

结果：

1
2
3
4
5
6
10.684
Married-civ-spouse
Married-civ-spouse    181
Never-married          18
Divorced                1
Name: count, dtype: int64

基于拉普拉斯机制，报告噪声最大值算法如下：

  • 对于每个$r\in\mathcal R$，计算噪声分数$u(x,r)+\operatorname{Lap}\left(\dfrac{\Delta u}{\epsilon}\right)$。
  • 输出噪声分数最大的元素$r\in\mathcal R$。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
import pandas,numpy
def score(data,option):
    return data.value_counts()[option]/1000
def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def report_noisy_max(x,R,u,sensitivity,epsilon):
    scores=[u(x,r) for r in R] #计算R中每个回复的分数
    noisy_scores=[laplace_mech(score,sensitivity,epsilon) for score in scores] #为每个分数增加噪声
    max_idx=numpy.argmax(noisy_scores) #找到最大分数对应的回复索引号
    return R[max_idx] #返回此索引号对应的回复
adult=pandas.read_csv("adult_with_pii.csv")
options=adult['Marital Status'].unique()
print(report_noisy_max(adult['Marital Status'],options,score,1,1))

r=[report_noisy_max(adult['Marital Status'],options,score,1,1) for i in range(200)]
print(pandas.Series(r).value_counts())

结果：

1
2
3
4
Married-civ-spouse
Married-civ-spouse    197
Never-married           3
Name: count, dtype: int64

稀疏向量技术

稀疏向量技术SVT可节省隐私消耗量，适用于在数据集上执行敏感度为$1$的问询流。稀疏向量技术只发布问询流中第一个通过测试的问询索引号，优势在于无论总共收到多少问询，消耗的总隐私消耗量时固定的。

稀疏向量技术最基础的实例是高于阈值算法。算法输入敏感度为$1$的问询流、数据集$\mathcal D$、阈值$T$和隐私参数$\epsilon$，算法满足$\epsilon$-差分隐私。下面场景基于选择求和问询的裁剪边界：获得多个不同裁剪边界并分别计算噪声裁剪边界，选择一个尽可能低且不会导致最终回复结果改变太大的一个裁剪边界。

1
2
3
4
5
6
7
def above_threshold(queries,df,T,epsilon):
    T_hat=T+numpy.random.laplace(loc=0,scale=2/epsilon) #噪声阈值  
    for idx,q in enumerate(queries):
        nu_i=numpy.random.laplace(loc=0,scale=4/epsilon)
        if q(df)+nu_i>=T_hat:
            return idx #第一个噪声问询回复大于噪声阈值的问询索引号
    return random.randint(0,len(queries)-1) #返回错误的索引号

高于阈值算法之所以满足差分隐私，是因为有可能返回错误的索引号，索引号对应的问询回复结果可能未超过给定阈值，索引号对应的问询有可能不是第一个回复结果超过阈值的问询。

对于上述选择求和问询的裁剪边界的场景，使用稀疏向量技术能获得更好的效果。

1
2
3
4
5
import pandas,numpy
def age_sum_query(df,b):
    return df['Age'].clip(lower=0,upper=b).sum()
adult=pandas.read_csv("adult_with_pii.csv")
print(age_sum_query(adult,30))

结果为：

1
879648

接下来定义一个求和差问询流，并基于稀疏向量技术，应用高于阈值确定b的最佳取值的问询索引号。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
def auto_avg(df,epsilon):
    def create_query(b):
        return lambda df:df.clip(lower=0,upper=b).sum()-df.clip(lower=0,upper=b+1).sum()
    bs=range(1,150000,5) #构造问询流
    queries=[create_query(b) for b in bs]
    epsilon_svt=epsilon/3 #使用1/3的隐私预算执行高于阈值 得到一个好的裁剪参数
    final_b=bs[above_threshold(queries,df,0,epsilon_svt)]
    epsilon_sum=epsilon/3 #分别使用1/3的隐私预算来获得噪声求和值与噪声计数值
    epsilon_count=epsilon/3
    noisy_sum=laplace_mech(df.clip(lower=0,upper=final_b).sum(),final_b,epsilon_sum)
    noisy_count=laplace_mech(len(df),1,epsilon_count)
    return noisy_sum/noisy_count
print(auto_avg(adult['Age'],1))
print(auto_avg(adult['Capital Gain'],1)) #可对不同尺度的数据使用相同的auto_avg函数

输出：

1
2
41.773740796468665
1113.791301048865

无论备选列表有多长，都能获得准确的结果，并消耗相同的隐私预算。该算法用了三次差分隐私机制，一次高于阈值算法，两次拉普拉斯机制，每个机制消耗$\dfrac{1}{3}\epsilon$的隐私预算，此算法满足$\epsilon$-差分隐私。

上述场景中只返回第一个超过阈值的问询索引号，有些场景中需要找到所有超过阈值的问询索引号。该任务可用稀疏算法挖成，但必须消耗更高的隐私预算。实现方法为：

  • 从问询流$q=\lbrace q_1,q_2,\cdots,q_k\rbrace$开始。
  • 在问询流$q$上执行高于阈值算法，得到第一个超过阈值的问询索引号$i$。
  • 用$q’=\lbrace q_{i+1},\cdots,q_k\rbrace$重启算法，转到第一步。
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
import pandas,numpy
def age_sum_query(df,b):
    return df['Age'].clip(lower=0,upper=b).sum()
def create_query(b):
    return lambda df:age_sum_query(df,b)-age_sum_query(df,b+1)
def above_threshold_fail_signal(queries,df,T,epsilon):
    T_hat=T+numpy.random.laplace(loc=0,scale = 2/epsilon)   
    for idx,q in enumerate(queries):
        nu_i=numpy.random.laplace(loc=0,scale = 4/epsilon)
        if q(df)+nu_i>=T_hat:
            return idx
    return None 
def sparse(queries,df,c,T,epsilon):
    idxs=[]
    pos=0
    epsilon_i=epsilon/c
    while pos<len(queries) and len(idxs)<c: #如果我们执行完问询流中的所有问询 或者我们找到了c个超过阈值的问询回复 则停止
        next_idx=above_threshold_fail_signal(queries[pos:],df,T,epsilon_i) #执行高于阈值 寻找下一个超过阈值的问询回复
        if next_idx==None: #如果高于阈值执行完了最后一个问询 则返回所有超过阈值的问询索引号
            return idxs
        pos=next_idx+pos #否则 使pos指向问询流中剩余的问询
        idxs.append(pos) #添加高于阈值找到的问询索引号
        pos=pos+1 #移动到问询流中的下一个问询
    return idxs
adult=pandas.read_csv("adult_with_pii.csv")
bs=range(1,150,5)
queries=[create_query(b) for b in bs]
epsilon=1
print(sparse(queries,adult,3,0,epsilon))

结果：

1
[15, 17, 18]

稀疏算法每次消耗$\dfrac{\epsilon}{c}$的隐私预算来调用高于阈值算法，满足$\epsilon$-差分隐私。

范围问询问的是“数据集中有多少行的值落在范围$(a,b)$中？”，是一种计数问询，敏感度为$1$。随机生成一组针对年龄列的范围问询，回复结果可能相差甚远，部分询问可能值匹配很少的数据行，计数值很小，另一部分可能匹配大量数据行，计数值很大。多数情况下，小计数值的差分隐私回复结果会很不准确，实际意义不大。要做的事了解那些问询结果事有价值的，并仅为这些偶价值的问询结果支付隐私预算。

此时确定一个阈值，并得到范围问询流中回复结果超过此阈值的问询索引号，认为这些问询是有价值的问询，随后用拉普拉斯机制得到这些有价值问询的差分隐私回复结果。

例如使用一般隐私预算来确定阈值为$10000$的前$c$个问询，另一半隐私预算用于获取这些问询的噪声回复结果。若高于阈值的问询数量远小于总问询数量，使用此方法可获得更准确的回复结果。

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
import pandas,numpy
def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def age_range_query(df,lower,upper):
    df1=df[df['Age']>lower]
    return len(df1[df1['Age']<upper])
def create_age_range_query():
    lower=numpy.random.randint(30,50)
    upper=numpy.random.randint(lower,70)
    return lambda df:age_range_query(df,lower,upper)
def above_threshold_fail_signal(queries,df,T,epsilon):
    T_hat=T+numpy.random.laplace(loc=0,scale = 2/epsilon)   
    for idx,q in enumerate(queries):
        nu_i=numpy.random.laplace(loc=0,scale = 4/epsilon)
        if q(df)+nu_i>=T_hat:
            return idx
    return None 
def sparse(queries,df,c,T,epsilon):
    idxs=[]
    pos=0
    epsilon_i=epsilon/c
    while pos<len(queries) and len(idxs)<c:
        next_idx=above_threshold_fail_signal(queries[pos:],df,T,epsilon_i)
        if next_idx==None:
            return idxs
        pos=next_idx+pos
        idxs.append(pos)
        pos=pos+1
    return idxs
def range_query_svt(queries,df,c,T,epsilon):
    sparse_epsilon=epsilon/2 #得到有价值的问询
    indices=sparse(queries,adult,c,T,sparse_epsilon)
    laplace_epsilon=epsilon/(2*c) #为每个有价值的问询执行拉普拉斯机制
    results=[laplace_mech(queries[i](df),1,laplace_epsilon) for i in indices]
    return results
adult=pandas.read_csv("adult_with_pii.csv")
range_queries=[create_age_range_query() for i in range(10)]
results=[q(adult) for q in range_queries]
print(results)
print(range_query_svt(range_queries,adult,5,10000,1))

结果为：

1
2
[0, 5798, 4748, 7870, 13678, 8474, 13678, 11064, 16800, 11558]
[13678.488637395474, 13684.047786305475, 11061.639560389842, 16810.671437590594, 11554.188968324295]

机器学习

一类监督学习问题：给定一组带标签的训练样本$\lbrace(x_1,y_1),(x_2,y_2),\cdots,(x_n,y_n)\rbrace$，其中$x_i$称为特征向量，$y_i$为标签，要训练一个模型$\theta$。该模型可预测没有在训练集中出现过的新特征向量对应的标签。一般每个$x_i$都是一个描述训练样本特征的实数向量，$y_i$是从预先定义好的类型集合中选取的，每个类型一般用一个整数来表示，要预先从全部样本中提取所有可能的类型，构成类型合集。二分类器的类型集合应包含两个类型，用$0/1$或$1/-1$表示。

训练模型时，从所有可用数据中选择一些来构造一组训练样本，留出一些作为测试样本。训练完模型应试验在非训练样本上的表现如何。若一个模型在未知新样本上表现很好，则称其泛化能力号。泛化能力不去的模型称其在训练数据上发生了过拟合。

给定一个包含$k$维特征向量$x_1,x_2,\cdots,x_k$的无标签样本，线性模型预测此样本的标签时，将计算$\displaystyle\sum_{i=1}^kw_ix_i+\beta$并用其符号作为预测标签，附属为$-1$，正数为$1$。$w_1,w_2,\cdots,w_k$称为模型的权重（或系数），$\beta$被称为偏差项（或截距），后续不考虑训练偏差项。

构建一个二分器的简单方法是使用逻辑回归，下面训练集包含$80\%$的样本，测试集包含剩余$20\%$的样本。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
import numpy
from sklearn.linear_model import LogisticRegression
X=numpy.load('adult_processed_x.npy')
y=numpy.load('adult_processed_y.npy')
training_size=int(X.shape[0]*0.8)
X_train=X[:training_size]
X_test=X[training_size:]
y_train=y[:training_size]
y_test=y[training_size:]
print(y_test.shape)
model=LogisticRegression().fit(X_train[:1000],y_train[:1000])
print(model.predict(X_test))
print(numpy.sum(model.predict(X_test)==y_test)/X_test.shape[0]) #准确率

def predict(xi,theta,bias=0): #手动构造线性模型
    label=numpy.sign(xi@theta+bias)
    return label
print(numpy.sum(predict(X_test,model.coef_[0],model.intercept_[0])==y_test)/X_test.shape[0]) #获得权重和偏差项

结果：

1
2
3
4
(9044,)
[-1. -1. -1. ... -1. -1. -1.]
0.8243034055727554
0.8243034055727554

损失函数衡量模型预测结果有多差，训练算法的目标是使损失函数达到最小值，损失值低的模型具有更好的预测能力。对于每个预测正确的样本，简单损失函数返回$0$，预测错误的返回$1$。二分类器常用对率损失函数，度量模型还有多远距离才能正确预测出标签。下面尝试对率损失函数，但权重全部简单设为$0$，所以样本全部预测错误，整个训练集平均损失值等于单个样本损失值。

1
2
3
4
5
6
def loss(theta,xi,yi):
    exponent=-yi*(xi.dot(theta))
    return numpy.log(1+numpy.exp(exponent))
theta=numpy.zeros(X_train.shape[1])
print(loss(theta,X_train[0],y_train[0]))
print(numpy.mean([loss(theta,x_i,y_i) for x_i,y_i in zip(X_train,y_train)]))

结果：

1
2
0.6931471805599453
0.6931471805599453

为了降低损失值，采用负梯度修改模型，这一过程称为梯度下降。

1
2
3
4
5
6
7
8
9
def gradient(theta,xi,yi): #梯度函数
    exponent=yi*(xi.dot(theta))
    return -(yi*xi)/(1+numpy.exp(exponent))
theta=theta-gradient(theta,X_train[0],y_train[0]) #本例只计算了训练集中第一个样本的梯度
print(theta)
print(y_train[0],predict(theta,X_train[0])) #更新后的模型具备分类样本的能力
def accuracy(theta): #度量准确性
    return numpy.sum(predict(X_test,theta)==y_test)/X_test.shape[0]
print(accuracy(theta))

结果：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
[ 0.          0.          0.          0.         -0.5         0.
  0.          0.         -0.5         0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.         -0.5         0.          0.          0.          0.
  0.          0.         -0.5         0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.         -0.5         0.          0.          0.
  0.          0.          0.          0.          0.          0.
 -0.5         0.         -0.5         0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.          0.
  0.          0.          0.          0.          0.         -0.5
  0.          0.         -0.25       -0.0606146  -0.21875     0.
  0.         -0.17676768]
-1.0 -1.0
0.7585139318885449

接下来把单样本梯度替换成所有训练样本的平均梯度，再定义一个可多次执行梯度下降的迭代算法，这就是常用的梯度下降算法。

1
2
3
4
5
6
7
8
9
10
11
def avg_grad(theta,X,y):
    grads=[gradient(theta,xi,yi) for xi,yi in zip(X,y)]
    return numpy.mean(grads,axis=0)
def gradient_descent(iterations):
    theta=numpy.zeros(X_train.shape[1])
    for i in range(iterations):
        theta=theta-avg_grad(theta,X_train,y_train)
    return theta
avg_grad(theta,X_train,y_train)
theta=gradient_descent(10) #可增大 100达到82%左右
print(accuracy(theta))

结果：

1
0.7787483414418399

接下来要设计一种算法为训练数据提供差分隐私保护，使最终训练得到的模型不会泄露与单个训练样本相关的任何信息。噪声梯度下降是在每轮模型更新前在梯度上增加噪声。

梯度是均值问询的结果，是每个样本梯度的均值。最好拆成一个求和问询和一个计数问询，这里计算每个样本梯度噪声和，除以噪声梯度值。此外需要限制每个样本梯度的敏感度，这里用类似采样—聚合框架的方法，裁剪梯度函数输出值，强制限定敏感度上限，称为梯度裁剪。

为了使输出向量的L2范数落在期望的范围内，通过缩放向量来做到裁剪向量，例如把向量中元素都除以向量的L2范数，则所得向量的L2范数为$1$。若将梯度表示为$\nabla(\theta;X,y)$，则梯度敏感度为$\left|\left|\operatorname{L2\_clip}(\nabla(\theta;X,y),b)-\operatorname{L2\_clip}(\nabla(\theta;X’,y),0)\right|\right|_2$，最差情况下$\left|\left|\operatorname{L2\_clip}(\nabla(\theta;X,y),b)\right|\right|_2=b,\operatorname{L2\_clip}(\nabla(\theta;X’,y),0)=\overrightarrow0$，两者L2范数差为$b$，做到了用裁剪参数$b$限定了梯度的L2敏感度上界。

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
def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def gaussian_mech_vec(v,sensitivity,epsilon,delta):
    return v+numpy.random.normal(loc=0,scale=sensitivity*numpy.sqrt(2*numpy.log(1.25/delta))/epsilon,size=len(v))
def L2_clip(v,b):
    norm = numpy.linalg.norm(v,ord=2)
    if norm>b:
        return b*(v/norm)
    else:
        return v
def gradient_sum(theta,X,y,b):
    gradients=[L2_clip(gradient(theta,x_i,y_i),b) for x_i,y_i in zip(X,y)]
    return numpy.sum(gradients,axis=0)
def noisy_gradient_descent(iterations,epsilon,delta):
    theta=numpy.zeros(X_train.shape[1])
    sensitivity=5.0
    noisy_count=laplace_mech(X_train.shape[0],1,epsilon)
    for i in range(iterations):
        grad_sum=gradient_sum(theta,X_train,y_train,sensitivity)
        noisy_grad_sum=gaussian_mech_vec(grad_sum,sensitivity,epsilon,delta)
        noisy_avg_grad=noisy_grad_sum/noisy_count
        theta=theta-noisy_avg_grad
    return theta
theta=noisy_gradient_descent(10,0.1,1e-5)
print(accuracy(theta))

结果：

1
0.7793011941618753

该算法每轮迭代满足$(\epsilon,\delta)$-差分隐私，且额外执行一次噪声计数问询得到满足$\epsilon$-差分隐私的噪声计数值。若执行$k$轮迭代，算法满足$((k+1)\epsilon,k\delta)$-差分隐私，可用高级组合性来分析总隐私消耗量。还可以将算法转为Rényi差分隐私或零集中差分隐私，得到更紧致的总隐私消耗量。

对率损失函数等梯度函数是Lipschitz连续的，这些梯度函数的全局敏感度是有界的。即如果$||x_i||_2\leqslant b$，则$||\nabla(\theta;x_i,y_i)||_2\leqslant b$。这样可通过裁剪训练样本，即梯度函数输入，来获得梯度函数的L2敏感度上界，不需要裁剪梯度函数输出。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def gradient_sum(theta,X,y,b):
    gradients=[gradient(theta,x_i,y_i) for x_i,y_i in zip(X,y)]
    return numpy.sum(gradients,axis=0)
def noisy_gradient_descent(iterations,epsilon,delta):
    theta=numpy.zeros(X_train.shape[1])
    sensitivity=5.0
    noisy_count=laplace_mech(X_train.shape[0],1,epsilon)
    clipped_X=[L2_clip(x_i,sensitivity) for x_i in X_train]
    for i in range(iterations):
        grad_sum=gradient_sum(theta,clipped_X,y_train,sensitivity)
        noisy_grad_sum=gaussian_mech_vec(grad_sum,sensitivity,epsilon,delta)
        noisy_avg_grad=noisy_grad_sum/noisy_count
        theta=theta-noisy_avg_grad
    return theta
theta=noisy_gradient_descent(10,0.1,1e-5)
print(accuracy(theta))

结果：

1
0.7826183104820875

更小的$\epsilon$会带来准确率更低的模型。更多的迭代次数意味更大的隐私消耗量，而在标准梯度下降算法中，更多迭代次数以为产出更好的模型。在差分隐私保护下，当总隐私预算保持不变时，更多迭代次数可能导致模型变得更糟糕，不得不使用更小的$\epsilon$支持更多轮迭代，这带来更大的噪声。

本地差分隐私

在中心模型中，原始敏感数据被汇总到单个数据集中，假定分析者是恶意的，单存在一个可信数据管理者，由他持有数据集并能正确执行分析者指定的差分隐私机制。在本地模型中，数据在离开数据主题控制前就已经满足差分隐私，如在将数据发送给数据管理者之前，用户就在自己设备上为自己数据添加噪声，数据管理者不需要可信，因为收集的是已经满足差分隐私的数据。

但在相同隐私预算下，对于相同问询，本地模型问询结果准确性通常比中心模型低几个数量级，意味着只有较少类型的问询适用于本地差分隐私。只有当数据量较大，参与者数量较多，差分隐私本地模型分析结果准确率才可满足实际要求。

随机应答是一种本地差分隐私机制，允许用户用错误的回复来应答调研中的敏感问题，引入的不确定性是差分隐私机制可提供隐私保护的根本原因，但算法输出结果仍释放足够信号帮助推断相关信息。在此过程中，数据主体按下述方法用“是”或“不是”回答一个问题：

  • 掷一枚硬币。
  • 若硬币正面向上，如实回答问题。
  • 若硬币反面向上，再掷一枚硬币。
  • 若第二枚硬币也是正面向上，回答“是”，否则回答“不是”。

该随机应答算法满足$\epsilon$-差分隐私，其中$\epsilon=\ln3$。

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
import numpy,pandas
def rand_resp_sales(response):
    truthful_response=response=='Sales'
    if numpy.random.randint(0,2)==0:
        return truthful_response
    else:
        return numpy.random.randint(0,2)==0
print(pandas.Series([rand_resp_sales('Sales') for i in range(200)]).value_counts())

adult=pandas.read_csv("adult_with_pii.csv")
responses=[rand_resp_sales(r) for r in adult['Occupation']]
print(pandas.Series(responses).value_counts())
print(len(adult[adult['Occupation']=='Sales'])) #“是”的数量

responses=[rand_resp_sales(r) for r in adult['Occupation']]
fake_yesses=len(responses)/4 #1/4的“是”来自随机结果 都是假的
num_yesses=numpy.sum([1 if r else 0 for r in responses]) #“是”总人数
true_yesses=num_yesses-fake_yesses #真实“是”的人数
rr_result = true_yesses*2
print(rr_result)

def laplace_mech(v,sensitivity,epsilon):
    return v+numpy.random.laplace(loc=0,scale=sensitivity/epsilon)
def pct_error(orig,priv):
    return numpy.abs(orig-priv)/orig*100.0
true_result=numpy.sum(adult['Occupation']=='Sales')
print(true_result)
print(pct_error(true_result,rr_result)) #错误率
print(pct_error(true_result,laplace_mech(true_result,1,1))) #中心模型拉普拉斯机制错误率

结果：

1
2
3
4
5
6
7
8
9
10
11
True     148
False     52
Name: count, dtype: int64
False    22735
True      9828
Name: count, dtype: int64
3650
3606.5
3650
1.1917808219178083
0.006456525121655132

一元编码方法解决本地差分隐私的直方图问询问题。首先定义应大于，即直方图包含的标签。然后下面三个函数实现了一元编码机制：

  • 编码函数encode：编码应答值。
  • 扰动函数perturb：扰动编码后的应答值。
  • 聚合函数aggregate：根据扰动应答值重构最终结果。

若应答域大小为$k$，将每个应答值编码为长度为$k$的比特向量，除应答者只为对应比特值为$1$以外，其他未知编码均为$0$，这种方法称为独热编码。

用扰动函数反转应答向量中各个比特值，满足差分隐私。翻转一个比特位的概率由$p$和$q$共同决定，即$\operatorname{Pr}(B’[i]=1)=\begin{cases}p&:B[i]=1\\\\q&:B[i]=0\end{cases}$。$\epsilon$可通过$p$和$q$计算得到$\epsilon=\ln\dfrac{p(1-q)}{(1-p)q}$。

聚合函数需要考虑每个标签的“假”应答数量，得到聚合结果$A[i]=\dfrac{\displaystyle\sum_j B_j’[i]-nq}{p-q}$。

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
import pandas,numpy
adult=pandas.read_csv("adult_with_pii.csv")
domain=adult['Occupation'].dropna().unique()
print(domain)

def encode(response):
    return [1 if d == response else 0 for d in domain]
print(encode('Sales'))

def perturb_bit(bit):
    p=.75
    q=.25
    sample=numpy.random.random()
    if bit==1:
        if sample<=p:
            return 1
        else:
            return 0
    elif bit==0:
        if sample<=q:
            return 1
        else: 
            return 0
def perturb(encoded_response):
    return [perturb_bit(b) for b in encoded_response]
print(perturb(encode('Sales')))

def unary_epsilon(p,q):
    return numpy.log((p*(1-q))/((1-p)*q))
print(unary_epsilon(.75,.25))

def aggregate(responses):
    p=.75
    q=.25
    sums=numpy.sum(responses,axis=0)
    n=len(responses)
    return [(v-n*q)/(p-q) for v in sums]
responses=[perturb(encode(r)) for r in adult['Occupation']]
counts=aggregate(responses)
print(list(zip(domain,counts)))

结果：

1
2
3
4
5
6
7
8
['Adm-clerical' 'Exec-managerial' 'Handlers-cleaners' 'Prof-specialty'
 'Other-service' 'Sales' 'Craft-repair' 'Transport-moving'
 'Farming-fishing' 'Machine-op-inspct' 'Tech-support' 'Protective-serv'
 'Baby' 'Armed-Forces' 'Priv-house-serv']
[0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[0, 1, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0]
2.1972245773362196
[('Adm-clerical', 3518.5), ('Exec-managerial', 3718.5), ('Handlers-cleaners', 1432.5), ('Prof-specialty', 4076.5), ('Other-service', 3512.5), ('Sales', 3630.5), ('Craft-repair', 4174.5), ('Transport-moving', 1630.5), ('Farming-fishing', 780.5), ('Machine-op-inspct', 1806.5), ('Tech-support', 1076.5), ('Protective-serv', 766.5), ('Baby', -217.5), ('Armed-Forces', -45.5), ('Priv-house-serv', 68.5)]

合成数据

合成数据生成算法的输入是一个原始数据集，输出是维度相同的合成数据集，期望合成数据集的数据与原始数据集的对应数据满足相同性质，且保留列之间的相关性。合成表示与原始数据维数不同，但可用于回答原始数据的问询，如直方图等。

为了满足差分隐私，在直方图每个计数值上单独增加拉普拉斯噪声即可，满足$\epsilon$-差分隐私。

对于拉普拉斯机制，计数值随着年龄范围增大而增大，更大分组意味着更强的信号，统计结果信噪比贬低，相对误差降低。对于合成表示，相对误差与问询范围无关，因为信号变强了，很多小分组的噪声也加在一起，噪声也变大了。

合成表示的下一步是合成数据，方法是将合成表示视为一个可用于估计原始数据潜在分布的概率分布函数，进而根据此概率分布采样，得到合成数据集。这里只考虑单列数据，其概率称为单维边际分布。先对直方图每个属性值计数，随后归一化计数结果，使所有计数结果和为$1$，得到概率值。

为了生成包含更多数据列的合成数据集，如考虑年龄和职业两列，可采用二维边缘分布。沿用之前的做法，在各个计数值上增加噪声后归一化处理，得到的计数结果视为概率值，每个计数结果对应的是年龄和职业一对属性值。

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
import pandas,numpy
def range_query(df,col,a,b):
    return len(df[(df[col]>=a)&(df[col]<b)])
def range_query_synth(syn_rep,a,b):
    total=0
    for i in range(a,b):
        total+=syn_rep[i]
    return total
adult=pandas.read_csv("adult_with_pii.csv")
bins=list(range(0, 100))
counts=[range_query(adult,'Age',b,b+1) for b in bins]
print(range_query_synth(counts,21,33))

def laplace_mech(v, sensitivity, epsilon):
    return v + numpy.random.laplace(loc=0, scale=sensitivity / epsilon)
epsilon = 1
dp_syn_rep = [laplace_mech(c, 1, epsilon) for c in counts]
print(range_query_synth(dp_syn_rep, 21, 33))

def pct_error(orig, priv):
    return numpy.abs(orig - priv)/orig * 100.0
true_answer = range_query(adult, 'Age', 30, 31)
print('合成表示误差率：{}'.format(pct_error(true_answer, range_query_synth(dp_syn_rep, 30, 31))))
print('拉普拉斯机制误差率：{}'.format(pct_error(true_answer, laplace_mech(true_answer, 1, epsilon))))
true_answer = range_query(adult, 'Age', 30, 71)
print('合成表示误差率：{}'.format(pct_error(true_answer, range_query_synth(dp_syn_rep, 30, 71))))
print('拉普拉斯机制误差率：{}'.format(pct_error(true_answer, laplace_mech(true_answer, 1, epsilon))))

dp_syn_rep_nn = numpy.clip(dp_syn_rep, 0, None)
syn_normalized = dp_syn_rep_nn / numpy.sum(dp_syn_rep_nn)
def gen_samples(n):
    return numpy.random.choice(bins, n, p=syn_normalized) #加权采样
syn_data = pandas.DataFrame(gen_samples(5), columns=['Age'])
print(syn_data)

print('年龄均值合成结果：{}'.format(numpy.mean(syn_data['Age'])))
print('年龄均值真实结果：{}'.format(numpy.mean(adult['Age'])))
print('误差率：{}'.format(pct_error(numpy.mean(syn_data['Age']), numpy.mean(adult['Age']))))
n = laplace_mech(len(adult), 1, 1.0)
syn_data = pandas.DataFrame(gen_samples(int(n)), columns=['Age'])
print('年龄范围问询合成结果：{}'.format(range_query(syn_data, 'Age', 20, 65)))
print('年龄范围问询真实结果：{}'.format(range_query(adult, 'Age', 20, 65)))
print('误差率：{}'.format(pct_error(range_query(adult, 'Age', 20, 65), range_query(syn_data, 'Age', 20, 65))))

ct=pandas.crosstab(adult['Age'],adult['Occupation'])
dp_ct=ct.applymap(lambda x:max(laplace_mech(x,1,1),0))
dp_vals=dp_ct.stack().reset_index().values.tolist()
probs=[p for _,_,p in dp_vals]
vals=[(a,b) for a,b,_ in dp_vals]
probs_norm=probs/numpy.sum(probs)
print(list(zip(vals,probs_norm))[0])
indices=range(0,len(vals))
n=laplace_mech(len(adult),1,1.0)
gen_indices=numpy.random.choice(indices,int(n),p=probs_norm)
syn_data=[vals[i] for i in gen_indices]
syn_df=pandas.DataFrame(syn_data,columns=['Age','Occupation'])
print(syn_df.head())
real_answer=range_query(adult,'Age',20,30)
syn_answer=range_query(syn_df,'Age',20,30)
print('误差率：',pct_error(real_answer,syn_answer))

结果：

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
6245
6243.090943873088
合成表示误差率：0.34286383388197766
拉普拉斯机制误差率：0.3164829375822671
合成表示误差率：0.013231792885249637
拉普拉斯机制误差率：0.002105545219908062
   Age
0   34
1   48
2   68
3   60
4   21
年龄均值合成结果：46.2
年龄均值真实结果：41.77250253355035
误差率：9.583327849458119
年龄范围问询合成结果：23575
年龄范围问询真实结果：23574
误差率：0.004241961482989734
((1, 'Adm-clerical'), 1.080324441269666e-05)
   Age       Occupation
0   45            Sales
1   29  Exec-managerial
2   38            Sales
3   45     Tech-support
4   15  Farming-fishing
误差率： 0.42024832855778416