2020DCIC智慧海洋建设算法赛学习03-特征工程

序：
特征工程往往是算法比赛中最至关重要的一环，一个好的特征工程能够让你的分数有大幅的提升，而如何做好特征工程、从哪些方面入手构建特征就需要经验积累和学习TOP选手的优秀方案。通过学习TOP选手开源代码的特征工程部分，我们可以发现，对于智慧海洋这样一个包含时序和空间信息的赛题，通常可以从以下几方面来构造特征。

0 基本预处理

在进行特征工程之前，我们需要对原始数据做一些基本的预处理。
部分原始训练数据是这样的：

通过上一篇博客的数据分析，我们知道数据中没有缺失值，因此不用做缺失值填充。对于这份数据，我们需要从以下两方面进行预处理：

• 异常值处理
  关于异常值处理方法我们在上一篇博客中已经讨论过了，因此这里不多做赘述。
• 数据类型预处理
  需要对一些列的数据类型进行转换和处理：一是将类别型的type列映射为数值；二是转换数据格式，将x、y坐标和速度特征的数据格式转换为浮点型，方向数据格式转换为整型；三是将time特征转换为datetime格式，并抽取出年、月、日、时等特征。当然，为了便于操作，也可以对列名进行修改。

label_dic1 = {'拖网': 0, '围网': 1, '刺网': 2}
label_dic2 = {0: '拖网', 1: '围网', 2: '刺网'}
name_dic = {'渔船ID': 'id', '速度': 'v', '方向': 'direct', 'type': 'label'}

# 修改列名
train_df.rename(columns=name_dic, inplace=True)

# 将标签映射为数值
train_df['label'] = train_df['label'].map(label_dic1)

# 转换数据格式
cols = ['x', 'y', 'v']
for col in cols:
    train_df[col] = train_df[col].astype('float')
train_df['direct'] = train_df['direct'].astype('int')

# 将时间数据转换成datetime格式，并抽取出相关的时间特征
train_df['time'] = pd.to_datetime(train_df['time'], format='%m%d %H:%M:%S')
train_df['date'] = train_df['time'].dt.date
train_df['hour'] = train_df['time'].dt.hour
train_df['month'] = train_df['time'].dt.month
train_df['weekday'] = train_df['time'].dt.weekday

train_df.head()

预处理后的训练数据如下图所示：

1 根据赛题背景构造特征

特征工程没有统一的模板，需要针对具体问题具体分析，深入挖掘赛题背景的相关知识，同时也包括一些常识性的知识，由此来构建一些与预测标签相关度高的特征。
例如，在这个赛题中，可以由此构造以下一些特征：

• 各坐标点与特定点之间的距离
• 划分白天和黑夜
• 划分季度
• 划分速度、方向等级
• 计算方向变化率和偏移角变化率
• 计算x轴和y轴上的速度
• 计算曲率

1.1 各坐标点与特定点之间的距离

以特定点(6165599, 5202660)为基准，计算各坐标点与基准点之间的距离。这里选取的特定点是复赛数据中的最值，在其他包含空间信息的类似场景中，我们在构建特征时也可以考虑以某个具有特殊意义的点作为基准，构造这样的距离特征。

train_df['x_abs_diff'] = (train_df['x'] - 6165599).abs()
train_df['y_abs_diff'] = (train_df['y'] - 5202660).abs()
# 计算与基准点之间的距离
train_df['base_dist'] = (train_df['x_abs_diff']**2 + train_df['y_abs_diff']**2) ** 0.5
# 删去无用的特征
del train_df['x_abs_diff'], train_df['y_abs_diff']

1.2 划分白天和黑夜

渔船并非总是处在作业状态，在晚上时渔船通常会停靠在港口，而我们的目标是要通过渔船的作业轨迹来判断该船是在进行什么类型的作业，因此我们可以构造一个特征来标记白天和黑夜。
取5-20时为白天，标记为1，其余时间为黑夜，标记为0：

train_df['day_night'] = 0
train_df.loc[(train_df['hour']>5) & (train_df['hour']<20), 'day_night'] = 1

1.3 划分季度

在不同季度鱼群的生长和活跃情况有所不同，因此对应的渔船生产作业方式也存在差异，因此我们可以构造一个特征来标记不同的季度：

train_df['quarter'] = 0
train_df.loc[(train_df['month'].isin([1, 2, 3])), 'quarter'] = 1
train_df.loc[(train_df['month'].isin([4, 5, 6])), 'quarter'] = 2
train_df.loc[(train_df['month'].isin([7, 8, 9])), 'quarter'] = 3
train_df.loc[(train_df['month'].isin([10, 11, 12])), 'quarter'] = 4

1.4 划分速度、方向等级

根据上一篇博客的数据分析，我们发现不同作业类型的渔船速度和方向分布存在较显著的区别，因此我们可以将速度和方向划分成多个等级，并统计每条渔船的速度和方向落在各个等级中的数目，由此区分渔船多数情况下处在哪些等级中。

• 将速度划分为0-7共八个等级，并统计每条渔船速度落在各个等级的数目

tmp_df = train_df.copy()
tmp_df.rename(columns={'id': 'ship', 'direct': 'd'}, inplace=True)

# 将速度划分等级
def vlevel(v):
    if v < 0.1:
        return 0
    elif v < 0.5:
        return 1
    elif v < 1:
        return 2
    elif v < 2.5:
        return 3
    elif v < 5:
        return 4
    elif v < 10:
        return 5
    elif v < 20:
        return 6
    else:
        return 7

# 统计每条渔船各个速度等级的数目
def get_vlevel_cnt(df):
    df['vlevel'] = df['v'].apply(lambda x: vlevel(x))
    tmp = df.groupby(['ship', 'vlevel'], as_index=False)['vlevel'].agg({'vlevel_cnt': 'count'})
    tmp = tmp.pivot(index='ship', columns='vlevel', values='vlevel_cnt')
    
    new_col_name = ['vlevel_' + str(col) for col in tmp.columns.tolist()]
    tmp.columns = new_col_name
    tmp = tmp.reset_index()
    
    return tmp

c1 = get_vlevel_cnt(tmp_df)
c1.head()

得到的渔船不同速度等级数目如下：

• 将方向区间0~360进行16等分，并统计每条渔船方向落在各个区间的数目

# 将方向分为16等分
def direct_level(df):
    df['d16'] = df['d'].apply(lambda x: int((x/22.5)+0.5)%16 if not np.isnan(x) else np.nan)
    return df

def get_direct_level_cnt(df):
    df = direct_level(df)
    tmp = df.groupby(['ship', 'd16'], as_index=False)['d16'].agg({'d16_cnt': 'count'})
    tmp = tmp.pivot(index='ship', columns='d16', values='d16_cnt')
    
    new_col_name = ['d16_' + str(col) for col in tmp.columns.tolist()]
    tmp.columns = new_col_name
    tmp = tmp.reset_index()
    
    return tmp

c2 = get_direct_level_cnt(tmp_df)
c2.head()

得到的部分渔船不同方向区间数目如下：

1.5 计算方向变化率和偏移角变化率

假设一条渔船的轨迹是1->2->3->4，如上图所示，其中蓝色实线表示渔船两位置间的轨迹线，红色虚线表示渔船在各个位置的方向，用${\theta }_1、{\theta }_2、{\theta }_3、{\theta }_4$表示方向角，${\alpha }_{12}、{\alpha }_{23}$表示两轨迹线之间的夹角，称为偏移角。
由此可以计算方向角变化率：
$$\frac{|{\theta }_1-{\theta }_2|}{timedif{f}_{12}}$$
偏移角变化率：
$$\frac{|{\alpha }_{12}-{\alpha }_{23}|}{timedif{f}_{23}}$$
其中$timedif{f}_{12}$表示位置1、2之间的时间差，$timedif{f}_{23}$表示位置2、3之间的时间差，单位为秒。
这里方向角就是方向特征，偏移角可以通过计算两条轨迹线与水平轴夹角的差得到。

import math
import time

# 计算每条渔船轨迹的方向变化率和偏移角变化率
def get_direct_change_rate(df):
    tmp = df.copy()
    tmp.sort_values(['ship', 'time'], ascending=True, inplace=True)
    
    # 获取下一个时间和x、y坐标
    tmp['time_next'] = tmp.groupby('ship')['time'].shift(-1)
    tmp['x_next'] = tmp.groupby('ship')['x'].shift(-1)
    tmp['y_next'] = tmp.groupby('ship')['y'].shift(-1)
    # 填充结尾的空值
    tmp['x_next'] = tmp['x_next'].fillna(method='ffill')
    tmp['y_next'] = tmp['y_next'].fillna(method='ffill')
    
    # 求轨迹线与水平轴之间的夹角
    tmp['change_angle'] = (tmp['y_next']-tmp['y']) / (tmp['x_next']-tmp['x'])
    tmp['change_angle'] = np.arctan(tmp['change_angle']) / math.pi * 180
    tmp['change_angle_next'] = tmp.groupby('ship')['change_angle'].shift(-1)
    
    # 获取时间差
    tmp['time_diff'] = np.abs(tmp['time_next'] - tmp['time'])
    tmp['time_diff'] = tmp['time_diff'].fillna(method='ffill')
    
    # 计算偏移角
    tmp['change_angle_diff'] = np.abs(tmp['change_angle'] - tmp['change_angle_next'])
    # 将偏移角控制在180度范围内
    tmp.loc[tmp['change_angle_diff']>180, 'change_angle_diff'] = 360 - tmp.loc[tmp['change_angle_diff']>180, 'change_angle_diff']
    # 求偏移角的变化率(度/秒)
    tmp['change_angle_v'] = tmp.apply(lambda x: x['change_angle_diff']/x['time_diff'].total_seconds(), axis=1)
    
    # 获取下一个方向角和方向角差
    tmp['d_next'] = tmp.groupby('ship')['d'].shift(-1)
    tmp['d_diff'] = np.abs(tmp['d_next'] - tmp['d'])
    # 将方向差控制在180度范围内
    tmp.loc[tmp['d_diff']>180, 'd_diff'] = 360 - tmp.loc[tmp['d_diff']>180, 'd_diff']
    # 求方向角的变化率(度/秒)
    tmp['d_v'] = tmp.apply(lambda x: x['d_diff']/x['time_diff'].total_seconds(), axis=1)
    
    # 获取偏移角变化率和方向变化率的最大值
    tmp1 = tmp[['ship', 'change_angle_v', 'd_v']]
    change_angle_v_max = tmp1.groupby('ship')['change_angle_v'].agg([('change_angle_v_max', 'max')])
    change_angle_v_max = change_angle_v_max.reset_index()
    d_v_max = tmp1.groupby('ship')['d_v'].agg([('d_v_max', 'max')])
    d_v_max = d_v_max.reset_index()
    
    # 将构造的两个特征合并成一张表格
    tmp = change_angle_v_max.merge(d_v_max, on='ship', how='left')
    
    return tmp

c5 = get_direct_change_rate(tmp_df)
c5.head()

构造的特征如下：

这里只取了方向变化率和偏移角变化率的最大值特征，实际可以保留方向变化率和偏移角变化率这两个特征。

1.6 计算x轴和y轴上的速度

尽管已经有了总的速度，我们还可以分别计算渔船在x轴和y轴上的速度：

tmp_df = train_df.copy()
tmp_df.rename(columns={'id': 'ship'}, inplace=True)

# 计算x和y轴上的速度
tmp_df.sort_values(['ship', 'time'], ascending=True, inplace=True)
tmp_df['y_next'] = tmp_df.groupby('ship')['y'].shift(-1)
tmp_df['x_next'] = tmp_df.groupby('ship')['x'].shift(-1)
tmp_df['y_next'] = tmp_df['y_next'].fillna(method='ffill')
tmp_df['x_next'] = tmp_df['x_next'].fillna(method='ffill')
tmp_df['time_next'] = tmp_df.groupby('ship')['time'].shift(-1)
tmp_df['time_diff'] = np.abs(tmp_df['time_next'] - tmp_df['time'])
tmp_df['v_y'] = tmp_df.apply(lambda x: (x['y_next']-x['y']) / x['time_diff'].total_seconds(), axis=1)
tmp_df['v_x'] = tmp_df.apply(lambda x: (x['x_next']-x['x']) / x['time_diff'].total_seconds(), axis=1)

1.7 计算曲率

我们也可以计算渔船的轨迹在每一个位置点的曲率：

# 计算曲率
tmp_df['y_prev'] = tmp_df.groupby('ship')['y'].shift(1)
tmp_df['x_prev'] = tmp_df.groupby('ship')['x'].shift(1)
tmp_df['y_prev'] = tmp_df['y_prev'].fillna(method='bfill')
tmp_df['x_prev'] = tmp_df['x_prev'].fillna(method='bfill')
tmp_df['dist_prev'] = ((tmp_df['x']-tmp_df['x_prev'])**2 + (tmp_df['y']-tmp_df['y_prev'])**2) ** 0.5
tmp_df['dist_next'] = ((tmp_df['x']-tmp_df['x_next'])**2 + (tmp_df['y']-tmp_df['y_next'])**2) ** 0.5
tmp_df['dist_prev_next'] = ((tmp_df['x_prev']-tmp_df['x_next'])**2 + (tmp_df['y_prev']-tmp_df['y_next'])**2) ** 0.5
tmp_df['curvature'] = (tmp_df['dist_prev']+tmp_df['dist_next']) / tmp_df['dist_prev_next']

2 根据空间信息构造分箱特征

分箱特征就是将特征的取值区间划分成多个小区间，也就是把大箱子分割成多个小箱子，然后把每条数据装进对应的小箱子里。例如在上一节中将速度和方向特征划分多个等级实际也是分箱特征。
构造分箱特征的方法除了上一节中用到的按区间划分和等分外，还可以按分位数分箱。分位数是指将一个变量按概率分布范围分为多个等份的数值点，例如常用的四分位数，其第一四分位数就是将所有数值从小到大排列后位于25%的数字。
我们也可以将速度按分位数分箱，例如按200分位数分箱，就是将所有速度值按从小到大排列，每隔1/200个数就是一个分位数，即一个分割点，由此将速度值划分到200个箱子中。这种按分位数分箱的方式可以用pd.qcut()实现。在分箱后对每个箱子进行编码，就得到了分箱特征。

# 对速度进行200分位数分箱，例如将前1/200分入一个桶
train_df['v_bin'] = pd.qcut(train_df['v'], 200, duplicates='drop') 
# 分箱后进行映射编码
train_df['v_bin'] = train_df['v_bin'].map(dict(zip(train_df['v_bin'].unique(), range(train_df['v_bin'].nunique()))))

对于空间信息，我们常用的分箱方式是等分，将x坐标和y坐标分别等分，就可以把整个空间划分为多个矩形网格。如下图所示，每个区域的编号是x坐标的分箱编号与y坐标的分箱编号的组合。

def traj_to_bin(traj, x_min, x_max, y_min, y_max, row_bins, col_bins):
    # row_bins = (x_max - x_min) / 700，即将x坐标每700划分为一个区域
    # col_bins = (y_max - y_min) / 3000，即将y坐标每3000划分为一个区域
    x_bins = np.linspace(x_min, x_max, endpoint=True, num=row_bins+1)  # 注意row_bins+1是包括两端点的划分点个数
    y_bins = np.linspace(y_min, y_max, endpoint=True, num=col_bins+1)
    
    # 确定每个x坐标属于哪个区域，对x坐标进行分箱编码
    traj.sort_values(by='x', inplace=True)
    x_res = np.zeros((len(traj), ))
    j = 0
    for i in range(1, row_bins+1):
        low, high = x_bins[i-1], x_bins[i]
        while j < len(traj):
            # low-0.001是为了保证数值稳定性，因为linspace得到的划分点精度可能与给出的坐标精度不同
            if (traj['x'].iloc[j] <= high) & (traj['x'].iloc[j] > low-0.001):
                x_res[j] = i
                j += 1
            else:
                break
    traj['x_grid'] = x_res
    traj['x_grid'] = traj['x_grid'].astype('int')
    traj['x_grid'] = traj['x_grid'].apply(str)
    
    # 确定每个y坐标属于哪个区域，对y坐标进行分箱编码
    traj.sort_values(by='y', inplace=True)
    y_res = np.zeros((len(traj), ))
    j = 0
    for i in range(1, col_bins+1):
        low, high = y_bins[i-1], y_bins[i]
        while j < len(traj):
            # low-0.001是为了保证数值稳定性，因为linspace得到的划分点精度可能与给出的坐标精度不同
            if (traj['y'].iloc[j] <= high) & (traj['y'].iloc[j] > low-0.001):
                y_res[j] = i
                j += 1
            else:
                break
    traj['y_grid'] = y_res
    traj['y_grid'] = traj['y_grid'].astype('int')
    traj['y_grid'] = traj['y_grid'].apply(str)
    
    # 确定每个x、y坐标对属于哪个区域，x坐标编码与y坐标编码组合成区域编码
    traj['x_y_grid'] = [i + '_' + j for i, j in zip(traj['x_grid'].values.tolist(), traj['y_grid'].values.tolist())]
    
    traj.sort_values(by='time', inplace=True)
    
    return traj

x_min, x_max = train_df['x'].min(), train_df['x'].max()
y_min, y_max = train_df['y'].min(), train_df['y'].max()
row_bins = math.ceil((x_max - x_min) / 700)  # 向上取整
col_bins = math.ceil((y_max - y_min) / 3000)

train_df = traj_to_bin(train_df, x_min, x_max, y_min, y_max, row_bins, col_bins)

对地图进行网格划分和编码时，也可以采用第一篇博客中提到的Geohash7编码。

3 构造统计特征

统计特征是特征工程中常用的一种构造特征，通常在各种场景下构造统计特征都能够在一定程度上提高模型的学习效果。
常用的统计特征包括：计数特征、最大值、最小值、均值、中位数、众数、标准差、偏斜度等。

• 计数特征

# 统计每个区域的数据条数
def get_grid_count(traj):
    grid_cnt_df = traj.groupby('x_y_grid').count().reset_index()
    grid_cnt_df = grid_cnt_df[['x_y_grid', 'x']]
    grid_cnt_df.rename({'x': 'grid_cnt'}, axis=1, inplace=True)
    
    return grid_cnt_df

# 统计每个区域的渔船个数
def get_grid_id_count(traj):
    grid_id_cnt_df = traj.groupby('x_y_grid')['id'].nunique().reset_index()
    grid_id_cnt_df = grid_id_cnt_df[['x_y_grid', 'id']]
    grid_id_cnt_df.rename({'id': 'grid_id_cnt'}, axis=1, inplace=True)
    
    return grid_id_cnt_df 

• 其他基本统计特征

raw_cols = train_df.columns.tolist()

# 获得序列的起始值
def start(x):
    try:
        return x[0]
    except:
        return None

# 获得序列的结束值
def end(x):
    try:
        return x[-1]
    except:
        return None

# 获得序列的众数    
def mode(x):
    try:
        return pd.Series(x).value_counts().index[0]
    except:
        return None

# 获取每条渔船的dist_prev_move_bin和v_bin序列
for f in ['dist_prev_move_bin', 'v_bin']:
    train_df[f + '_seq'] = train_df['id'].map(train_df.groupby('id')[f].agg(lambda x: ','.join(x.astype(str))))
    
# 对序列构造一些基本的统计特征，其中np.ptp用于计算最大最小值之间的差值
tmp_df = train_df.groupby('id').agg({
    'id': ['count'], 'x_bin1': [mode], 'y_bin1': [mode], 'x_bin2': [mode], 'y_bin2': [mode], 'x_y_bin1': [mode],
    'x': ['mean', 'max', 'min', 'std', np.ptp, start, end], 'y': ['mean', 'max', 'min', 'std', np.ptp, start, end],
    'v': ['mean', 'max', 'min', 'std', np.ptp], 'direct': ['mean'], 'x_bin1_cnt': ['mean', 'max', 'min'], 'y_bin1_cnt': ['mean', 'max', 'min'],
    'x_bin2_cnt': ['mean', 'max', 'min'], 'y_bin2_cnt': ['mean', 'max', 'min'], 'x_y_bin1_cnt': ['mean', 'max', 'min'],
    'dist_prev_move': ['mean', 'max', 'min', 'std', 'sum'], 'y_xbin1_ymax_diff': ['mean', 'min'], 'y_xbin1_ymin_diff': ['mean', 'min'],
    'x_ybin1_xmax_diff': ['mean', 'min'], 'x_ybin1_xmin_diff': ['mean', 'min']
}).reset_index()

tmp_df.columns = ['_'.join(col).strip()for col in tmp_df.columns]
tmp_df.rename(columns={'id_': 'id'}, inplace=True)
cols = [f for f in tmp_df.keys() if f !='id']

我们可以对白天、黑夜、速度为零和速度非零的数据分别构造统计特征：

def group_feature(df, key, target, aggs, flag):
    agg_dic = []
    for ag in aggs:
        agg_dic.append(('{}_{}_{}'.format(target, ag, flag), ag))
    t = df.groupby(key)[target].agg(agg_dic).reset_index()    
    return t

def extract_feature(df, train, flag):
    if (flag=='on_night') or (flag=='on_day'):
        t = group_feature(df, 'ship', 'speed', ['max', 'median', 'mean', 'std', 'skew'], flag)
        train = pd.merge(train, t, on='ship', how='left')
    
    if flag == '0':
        t = group_feature(df, 'ship', 'direct', ['max', 'median', 'mean', 'std', 'skew'], flag)
        train = pd.merge(train, t, on='ship', how='left')
    elif flag == '1':
        t = group_feature(df, 'ship', 'speed', ['max', 'median', 'mean', 'std', 'skew'], flag)
        train = pd.merge(train, t, on='ship', how='left')
        t = group_feature(df, 'ship', 'direct', ['max', 'median', 'mean', 'std', 'skew'], flag)
        train = pd.merge(train, t, on='ship', how='left')
        speed_nunique = df.groupby('ship')['speed'].nunique().to_dict()
        train['speed_nunique_{}'.format(flag)] = train['ship'].map(speed_nunique)
        direct_nunique = df.groupby('ship')['direct'].nunique().to_dict()
        train['direct_nunique_{}'.format(flag)] = train['ship'].map(direct_nunique)
        
    t = group_feature(df, 'ship', 'x', ['max', 'min', 'median', 'mean', 'std', 'skew'], flag)
    train = pd.merge(train, t, on='ship', how='left')
    t = group_feature(df, 'ship', 'y', ['max', 'min', 'median', 'mean', 'std', 'skew'], flag)
    train = pd.merge(train, t, on='ship', how='left')
    t = group_feature(df, 'ship', 'base_dist', ['max', 'min', 'mean', 'std', 'skew'], flag)
    train = pd.merge(train, t, on='ship', how='left')
    
    train['x_max_x_min_{}'.format(flag)] = train['x_max_{}'.format(flag)] - train['x_min_{}'.format(flag)]
    train['y_max_y_min_{}'.format(flag)] = train['y_max_{}'.format(flag)] - train['y_min_{}'.format(flag)]
    train['y_max_x_min_{}'.format(flag)] = train['y_max_{}'.format(flag)] - train['x_min_{}'.format(flag)]
    train['x_max_y_min_{}'.format(flag)] = train['x_max_{}'.format(flag)] - train['y_min_{}'.format(flag)]
    # 计算斜率
    train['slope_{}'.format(flag)] = train['y_max_y_min_{}'.format(flag)] / np.where(train['x_max_x_min_{}'.format(flag)]==0, 0.001, train['x_max_x_min_{}'.format(flag)])
    # 计算面积
    train['area_{}'.format(flag)] = train['x_max_x_min_{}'.format(flag)] * train['y_max_y_min_{}'.format(flag)]
    
    mode_hour = df.groupby('ship')['hour'].agg(lambda x: x.value_counts().index[0]).to_dict()
    train['mode_hour_{}'.format(flag)] = train['ship'].map(mode_hour) 
    train['median_slope_{}'.format(flag)] = train['y_median_{}'.format(flag)] / np.where(train['x_median_{}'.format(flag)]==0, 0.001, train['x_median_{}'.format(flag)])
    
    return train

data = train_df.copy()
data.rename(columns={'id': 'ship', 'v': 'speed'}, inplace=True)
data_label = data.drop_duplicates(['ship'], keep='first')

# 对速度为零和速度非零的数据分别构造统计特征
data1 = data[data['speed']==0]
data2 = data[data['speed']!=0]
data_label = extract_feature(data1, data_label, '0')
data_label = extract_feature(data2, data_label, '1')

# 对白天和黑夜的数据分别构造统计特征
data1 = data[data['day_night']==0]
data2 = data[data['day_night']==1]
data_label = extract_feature(data1, data_label, 'on_night')
data_label = extract_feature(data2, data_label, 'on_day')

data_label.rename(columns={'ship': 'id', 'speed': 'v'}, inplace=True)

• 一些不常使用的统计特征

# 计算变异系数/离散系数
def coefficient_of_variation(x):
    x = x.values
    if np.mean(x) == 0:
        return 0
    return np.std(x) / np.mean(x)

# 获取第二大的数
def max2(x):
    x = list(x.values)
    x.sort(reverse=True)
    return x[1]

# 获取第三大的数
def max3(x):
    x = list(x.values)
    x.sort(reverse=True)
    return x[2]

# 统计绝对值均值
def diff_abs_mean(x):
    return np.mean(np.abs(np.diff(x)))

4 根据时序信息构造Embedding特征

在这个赛题中，每条渔船的轨迹是时序数据，它包括渔船的速度序列、方向序列和坐标序列，这里的坐标序列采用之前构造的区域编码序列。如何利用这些时序数据？一种常用的做法就是构造Embedding特征，Embedding能够将长短不一的序列映射成为一个等长的向量供模型进行学习。
构造Embedding特征通常采用两种方法：Word2vec和NMF。

4.1 采用Word2vec构造Embedding特征

Word2vec是自然语言处理中常用的一种算法，顾名思义，它的目标就是要将输入的语句映射成一个向量。
我们先定义一个用Word2vec构造词向量的函数：

# 构造渔船轨迹序列的词向量
def traj_cbow_embedding(traj_data_corpus=None, embedding_size=70, iters=40, min_count=3, window_size=25, seed=9012, num_runs=5, word_feat='x_y_grid'):
    ship_id = traj_data_corpus['id'].unique()
    sentences, embedding_df_list, embedding_model_list = [], [], []
    for i in ship_id:
        traj = traj_data_corpus[traj_data_corpus['id']==i]
        sentences.append(traj[word_feat].values.tolist())
        
    print('\n@Start CBOW word embedding at {}'.format(datetime.now()))
    print('-------------------------------------------')
    for i in tqdm(range(num_runs)):
        model = Word2Vec(sentences=sentences, vector_size=embedding_size, min_count=min_count, workers=mp.cpu_count(), window=window_size, seed=seed, epochs=iters, sg=0)
        embedding_vec = []
        for ind, seq in enumerate(sentences):
            seq_vec, word_count = 0, 0
            for word in seq:
                if word not in model.wv:
                    continue
                else:
                    seq_vec += model.wv[word]
                    word_count += 1
            if word_count == 0:
                embedding_vec.append(embedding_size * [0])
            else:
                embedding_vec.append(seq_vec / word_count)
        embedding_vec = np.array(embedding_vec)
        embedding_cbow_df = pd.DataFrame(embedding_vec, columns=['embedding_cbow_{}_{}'.format(word_feat, i) for i in range(embedding_size)])
        embedding_cbow_df['id'] = ship_id
        embedding_df_list.append(embedding_cbow_df)
        embedding_model_list.append(model)
    print('-------------------------------------------')
    print('\n@End CBOW word embedding at {}'.format(datetime.now()))
    return embedding_df_list, embedding_model_list

有了以上函数，我们就可以对速度序列、方向序列、坐标序列构造词向量，也可以构造速度和方向的组合特征序列的词向量：

# 构造渔船轨迹的坐标序列的词向量
embedding_size = 70
iters = 70
min_count = 3
window_size = 25
num_runs = 1
seed = 9012
word_feat = 'x_y_grid'

df_list, model_list = traj_cbow_embedding(train_df, embedding_size, iters, min_count, window_size, seed, num_runs, word_feat)
train_embedding_df_list = [df.reset_index(drop=True) for df in df_list]
feat = train_embedding_df_list[0]
feat = pd.DataFrame(feat)

ship_id = train_df['id'].unique()
total_embedding = pd.DataFrame(ship_id, columns=['id'])
traj_data = train_df[['v', 'direct', 'id']].rename(columns={'v': 'speed'})

# 构造速度、方向以及速度+方向组合特征序列
traj_data_corpus = []
traj_data['speed_str'] = traj_data['speed'].apply(lambda x: str(int(x*100)))
traj_data['direct_str'] = traj_data['direct'].apply(str)
traj_data['speed_direct_str'] = traj_data['speed_str'] + '_' + traj_data['direct_str']
traj_data_corpus = traj_data[['id', 'speed_str', 'direct_str', 'speed_direct_str']]

# 构造速度序列的词向量
print('\n@Round 2 speed embedding:')
df_list, model_list = traj_cbow_embedding(traj_data_corpus, embedding_size=10, iters=40, window_size=25, seed=9102, num_runs=1, word_feat='speed_str')
speed_embedding = df_list[0].reset_index(drop=True)
total_embedding = pd.merge(total_embedding, speed_embedding, on='id', how='left')

# 构造方向序列的词向量
print('\n@Round 3 direct embedding:')
df_list, model_list = traj_cbow_embedding(traj_data_corpus, embedding_size=10, iters=40, window_size=25, seed=9102, num_runs=1, word_feat='direct_str')
direct_embedding = df_list[0].reset_index(drop=True)
total_embedding = pd.merge(total_embedding, direct_embedding, on='id', how='left')

# 构造速度+方向的组合特征序列的词向量
print('\n@Round 4 speed—direct embedding:')
df_list, model_list = traj_cbow_embedding(traj_data_corpus, embedding_size=12, iters=70, window_size=25, seed=9102, num_runs=1, word_feat='speed_direct_str')
speed_direct_embedding = df_list[0].reset_index(drop=True)
total_embedding = pd.merge(total_embedding, speed_direct_embedding, on='id', how='left')

4.2 采用NMF算法构造Embedding特征

NMF，即非负矩阵分解，NMF算法对任一个给定的非负矩阵V(维度为MxN)，能够将其分解为非负矩阵W(维度为MxK)和非负矩阵H(维度为KxN)，即V=W*H，其中N为数据条数，M为每条数据的特征数，K为设定的主题个数，也就是说，对于一条长度为M的输入数据，NMF算法能够得到长度为K的输出向量，即将长度为M的序列映射为一个长度为K的向量。
需要注意的是，NMF算法要求的输入是等长的数值型序列，因此我们需要先将所有渔船的轨迹序列拼接成一整个文本，然后用TF-Idf算法将文本处理成一个MxN的矩阵，最后再用NMF算法得到Embedding特征。
我们先定义一个NMF类：

# 定义一个NMF类
class nmf_list(object):
    def __init__(self, data, by_name, to_list, nmf_n, top_n):
        self.data = data
        self.by_name = by_name
        self.to_list = to_list
        self.nmf_n = nmf_n
        self.top_n = top_n
        
    def run(self, tf_n):
        df_all = self.data.groupby(self.by_name)[self.to_list].apply(lambda x: '|'.join(x)).reset_index()
        self.data = df_all.copy()
        
        print('Build Word Frequency------------')
        # 统计词频
        def word_frec(x):
            word_dict = []
            x = x.split('|')
            docs = []
            for doc in x:
                doc = doc.split()
                docs.append(doc)
                word_dict.extend(doc)
            word_dict = Counter(word_dict)
            new_word_dict = {}
            for key, value in word_dict.items():
                new_word_dict[key] = [value, 0]
            del word_dict
            del x
            for doc in docs:
                doc = Counter(doc)
                for word in doc.keys():
                    new_word_dict[word][1] += 1
            return new_word_dict
        self.data['word_frec'] = self.data[self.to_list].apply(word_frec)
        
        print('Build Top_' + str(self.top_n) + '------------')
        # 获取前100个高频词列表
        def top_100(word_dict):
            return sorted(word_dict.items(), key = lambda x: (x[1][1], x[1][0]), reverse=True)[:self.top_n]
        self.data['top_' + str(self.top_n)] = self.data['word_frec'].apply(top_100)
        # 获取前100个高频词
        def top_100_word(word_list):
            words = []
            for i in word_list:
                i = list(i)
                words.append(i[0])
            return words
        self.data['top_' + str(self.top_n) + '_word'] = self.data['top_' + str(self.top_n)].apply(top_100_word)
        print(self.data.shape)
        
        # 统计前100个高频词的词频
        word_list = []
        for i in self.data['top_' + str(self.top_n) + '_word'].values:
            word_list.extend(i)
        word_list = Counter(word_list)
        word_list = sorted(word_list.items(), key=lambda x: x[1], reverse=True)
        user_frec = []
        for i in word_list:
            i = list(i)
            user_frec.append(i[1] / self.data[self.by_name].nunique())
        # 去掉停止词，即用户频率超过0.5的词，这是为了发掘每个用户的特殊性
        stop_words = []
        for i, j in zip(word_list, user_frec):
            if j > 0.5:
                i = list(i)
                stop_words.append(i[0])
        
        print('Start Title Feature------------')
        # 将去掉停止词的序列拼接成一个文本
        self.data['title_feature'] = self.data[self.to_list].apply(lambda x: x.split('|'))
        self.data['title_feature'] = self.data['title_feature'].apply(lambda line: [w for w in line if w not in stop_words])
        self.data['title_feature'] = self.data['title_feature'].apply(lambda x: ' '.join(x))
        
        print('Start NMF------------')
        # 采用TF-Idf算法处理文本，再用NMF算法提取文本主题
        tfidf_vectorizer = TfidfVectorizer(ngram_range=(tf_n, tf_n))
        tfidf = tfidf_vectorizer.fit_transform(self.data['title_feature'].values)
        text_nmf = NMF(n_components=self.nmf_n).fit_transform(tfidf)
        
        # 整理并输出
        name = [str(tf_n) + self.to_list + '_' + str(x) for x in range(1, self.nmf_n+1)]
        tag_list = pd.DataFrame(text_nmf)
        print(tag_list.shape)
        tag_list.columns = name
        tag_list[self.by_name] = self.data[self.by_name]
        column_name = [self.by_name] + name
        tag_list = tag_list[column_name]
        return tag_list

用NMF类构造速度、x和y序列的Embedding特征，这里我们将K设为8：

tmp_df = train_df.copy()
tmp_df.rename(columns={'v': 'speed', 'id': 'ship'}, inplace=True)
for j in range(1, 4):
    print('------------ {} ------------'.format(j))
    for i in ['speed', 'x', 'y']:
        tmp_df[i + '_str'] = tmp_df[i].astype(str)
        nmf = nmf_list(tmp_df, 'ship', i + '_str', 8, 100)
        nmf_a = nmf.run(j)
        nmf_a.rename(columns={'ship': 'id'}, inplace=True)
        data_label = data_label.merge(nmf_a, on='id', how='left')

得到Embedding特征如下图所示：

5 总结

以上就是本题中构造特征的一些方法，这些方法没有必要全部使用，在模型学习中也不一定都有效果，可以根据模型的学习结果进行筛选。同时，这些方法可以为其他题目进行特征工程提供一些可行的思路。