68 분 소요

데이터 분석 예제

드디어 마지막 장이다. 여기서는 실제 데이터셋을 살펴본다. 지금까지 책에서 배운 기술을 이용해서 데이터에서 의미 있는 정보를 추출해 보도록 하자. 여기서 설명하는 기술은 여러분의 데이터셋을 포함하여 모든 데이터셋에 적용할 수 있을 것이다. 이 장에는 이 책에서 배웠던 도구들을 실습해볼 수 있는 예제 데이터셋들을 모아두었다.

책에서 사용한 예제 데이터셋은 이 책의 깃허브 저장소에서 다운로드할 수 있다.

  • 깃허브 저장소

https://github.com/wesm/pydata-book

Bit.ly의 1.USA.gov 데이터

2011년 URL 축약 서비스인 Bit.ly는 미국 정부 웹사이트인 USA.gov와 제휴하여 .gov나 .mil로 끝나는 URL을 축약한 사용자들에 대한 익명 정보를 제공했었다. 2011년에는 실시간 피드뿐 아니라 매 시간마다 스냅샷을 텍스트 파일로 내려받을 수 있었다. 저자가 이 책을 쓰는 시기에 해당 서비스는 더 이상 존재하지 않지만 그 데이터 파일 중 하나를 살펴보자.

매 시간별 스냅샷 파일의 각 로우는 웹 데이터 형식으로 흔히 사용되는 JSONJavaScript Object Notation이다. 스냅샷 파일의 첫 줄을 열어보면 다음과 비슷한 내용을 확인할 수 있다.

path = 'datasets/bitly_usagov/example.txt'

open(path).readline()
'{ "a": "Mozilla\\/5.0 (Windows NT 6.1; WOW64) AppleWebKit\\/535.11 (KHTML, like Gecko) Chrome\\/17.0.963.78 Safari\\/535.11", "c": "US", "nk": 1, "tz": "America\\/New_York", "gr": "MA", "g": "A6qOVH", "h": "wfLQtf", "l": "orofrog", "al": "en-US,en;q=0.8", "hh": "1.usa.gov", "r": "http:\\/\\/www.facebook.com\\/l\\/7AQEFzjSi\\/1.usa.gov\\/wfLQtf", "u": "http:\\/\\/www.ncbi.nlm.nih.gov\\/pubmed\\/22415991", "t": 1331923247, "hc": 1331822918, "cy": "Danvers", "ll": [ 42.576698, -70.954903 ] }\n'

파이썬에는 JSON 문자열을 파이썬 사전 객체로 바꿔주는 다양한 내장 모듈과 서드파티 모듈이 있다. 여기서는 json 모듈의 loads() 함수를 이용해서 내려받은 샘플 파일을 한 줄씩 읽는다.

import json

path = 'datasets/bitly_usagov/example.txt'
records = [json.loads(line) for line in open(path, 'r', encoding='UTF-8')]

결과를 담고 있는 records 객체는 파이썬 사전의 리스트다.

records[0]
{'a': 'Mozilla/5.0 (Windows NT 6.1; WOW64) AppleWebKit/535.11 (KHTML, like Gecko) Chrome/17.0.963.78 Safari/535.11',
 'c': 'US',
 'nk': 1,
 'tz': 'America/New_York',
 'gr': 'MA',
 'g': 'A6qOVH',
 'h': 'wfLQtf',
 'l': 'orofrog',
 'al': 'en-US,en;q=0.8',
 'hh': '1.usa.gov',
 'r': 'http://www.facebook.com/l/7AQEFzjSi/1.usa.gov/wfLQtf',
 'u': 'http://www.ncbi.nlm.nih.gov/pubmed/22415991',
 't': 1331923247,
 'hc': 1331822918,
 'cy': 'Danvers',
 'll': [42.576698, -70.954903]}

순수 파이썬으로 표준시간대 세어보기

이 데이터에서 가장 빈도가 높은 표준시간대(tz 필드)를 구한다고 가정하자. 다양한 방법이 있지만 먼저 리스트 표기법을 사용해서 표준시간대의 목록을 가져오자.

time_zones = [rec['tz'] for rec in records]

하지만 records의 아이템이 모두 표준시간대 필드를 가지고 있는 건 아니라는 게 드러났다! 이 문제는 if ‘tz’ in rec을 리스트 표기법 뒤에 추가해서 tz 필드가 있는지 검사하면 쉽게 해결할 수 있다.

time_zones = [rec['tz'] for rec in records if 'tz' in rec]

time_zones[:10]
['America/New_York',
 'America/Denver',
 'America/New_York',
 'America/Sao_Paulo',
 'America/New_York',
 'America/New_York',
 'Europe/Warsaw',
 '',
 '',
 '']

상위 10개의 표준시간대를 보면 그중 몇 개는 비어 있어서 뭔지 알 수 없다. 비어 있는 필드를 제거할 수도 있지만 일단은 그냥 두고 표준시간대를 세어보자.

def get_counts(sequence):
    counts = {}
    for x in sequence:
        if x in counts:
            counts[x] += 1
        else:
            counts[x] = 1
    return counts

파이썬 표준 라이브러리에 익숙하다면 다음처럼 좀 더 간단하게 작성할 수도 있다.

from collections import defaultdict

def get_counts(sequence):
    counts = defaultdict(int) # 값이 0으로 초기화된다.
    for x in sequence:
        counts[x] += 1
    return counts

재사용이 쉽도록 이 로직을 함수로 만들고 이 함수에 time_zones 리스트를 넘겨서 사용하자.

counts = get_counts(time_zones)
counts['America/New_York']
1251
len(time_zones)
3440

가장 많이 등장하는 상위 10개의 표준시간대를 알고 싶다면 좀 더 세련된 방법으로 사전을 사용하면 된다.

def top_counts(count_dict, n=10):
    value_key_pairs = [(count, tz) for tz, count in count_dict.items()]
    value_key_pairs.sort()
    return value_key_pairs[-n:]

이제 상위 10개의 표준시간대를 구했다.

top_counts(counts)
[(33, 'America/Sao_Paulo'),
 (35, 'Europe/Madrid'),
 (36, 'Pacific/Honolulu'),
 (37, 'Asia/Tokyo'),
 (74, 'Europe/London'),
 (191, 'America/Denver'),
 (382, 'America/Los_Angeles'),
 (400, 'America/Chicago'),
 (521, ''),
 (1251, 'America/New_York')]

파이썬 표준 라이브러리의 collections.Counter 클래스를 이용하면 지금까지 했던 작업을 훨씬 쉽게 할 수 있다.

from collections import Counter

counts = Counter(time_zones)

counts.most_common(10)
[('America/New_York', 1251),
 ('', 521),
 ('America/Chicago', 400),
 ('America/Los_Angeles', 382),
 ('America/Denver', 191),
 ('Europe/London', 74),
 ('Asia/Tokyo', 37),
 ('Pacific/Honolulu', 36),
 ('Europe/Madrid', 35),
 ('America/Sao_Paulo', 33)]

pandas로 표준시간대 세어보기

records를 가지고 DataFrame을 만드는 방법은 아주 쉽다. 그냥 레코드가 담긴 리스트를 pandas.DataFrame으로 넘기면 된다.

import pandas as pd

frame = pd.DataFrame(records)

frame.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 3560 entries, 0 to 3559
Data columns (total 18 columns):
 #   Column       Non-Null Count  Dtype  
---  ------       --------------  -----  
 0   a            3440 non-null   object 
 1   c            2919 non-null   object 
 2   nk           3440 non-null   float64
 3   tz           3440 non-null   object 
 4   gr           2919 non-null   object 
 5   g            3440 non-null   object 
 6   h            3440 non-null   object 
 7   l            3440 non-null   object 
 8   al           3094 non-null   object 
 9   hh           3440 non-null   object 
 10  r            3440 non-null   object 
 11  u            3440 non-null   object 
 12  t            3440 non-null   float64
 13  hc           3440 non-null   float64
 14  cy           2919 non-null   object 
 15  ll           2919 non-null   object 
 16  _heartbeat_  120 non-null    float64
 17  kw           93 non-null     object 
dtypes: float64(4), object(14)
memory usage: 500.8+ KB
frame['tz'][:10]
0     America/New_York
1       America/Denver
2     America/New_York
3    America/Sao_Paulo
4     America/New_York
5     America/New_York
6        Europe/Warsaw
7                     
8                     
9                     
Name: tz, dtype: object

frame의 출력 결과는 거대한 DataFrame 객체의 요약 정보다. frame[‘tz’]에서 반환되는 Series 객체에는 value_counts 메서드를 이용해서 시간대를 세어볼 수 있다.

tz_counts = frame['tz'].value_counts()

tz_counts[:10]
America/New_York       1251
                        521
America/Chicago         400
America/Los_Angeles     382
America/Denver          191
Europe/London            74
Asia/Tokyo               37
Pacific/Honolulu         36
Europe/Madrid            35
America/Sao_Paulo        33
Name: tz, dtype: int64

matplotlib 라이브러리로 이 데이터를 그래프로 그릴 수 있다. 그전에 records에서 비어 있는 표준시간대를 다른 이름으로 바꿔보자. fillna() 함수로 빠진 값을 대체하고, 불리언 배열 색인을 이용해서 비어 있는 값을 대체할 수 있다.

clean_tz = frame['tz'].fillna('Missing')

clean_tz[clean_tz == ''] = 'Unknown'

tz_counts = clean_tz.value_counts()

tz_counts[:10]
America/New_York       1251
Unknown                 521
America/Chicago         400
America/Los_Angeles     382
America/Denver          191
Missing                 120
Europe/London            74
Asia/Tokyo               37
Pacific/Honolulu         36
Europe/Madrid            35
Name: tz, dtype: int64

여기서는 seaborn 패키지를 이용해서 수평막대그래프를 그려보자(그림 14-1).

# matplotlib 한글 폰트 오류 문제 해결
from matplotlib import font_manager, rc
font_path = "./malgun.ttf"    # 폰트 파일 위치
font_name = font_manager.FontProperties(fname=font_path).get_name()
rc('font', family=font_name)
import matplotlib.pyplot as plt
import seaborn as sns

subset = tz_counts[:10]

fig = plt.figure(figsize=(15, 10))
ax = fig.add_subplot(1, 1, 1)

sns.barplot(x=subset.values, y=subset.index, ax=ax)

plt.title('그림 14-1 1.usa.gov 예제 데이터에서 가장 많이 나타난 시간대', size=20, loc='left')
plt.xticks(fontsize=15)
plt.yticks(fontsize=15)

plt.show()

a 필드에는 URL 단축을 실행하는 브라우저, 단말기, 애플리케이션에 대한 정보(User Agent 문자열)가 들어 있다.

frame['a'][1]
'GoogleMaps/RochesterNY'
frame['a'][50]
'Mozilla/5.0 (Windows NT 5.1; rv:10.0.2) Gecko/20100101 Firefox/10.0.2'
frame['a'][51][:50] # 긴 문자열
'Mozilla/5.0 (Linux; U; Android 2.2.2; en-us; LG-P9'

‘agent’라고 하는 흥미로운 문자열 정보를 분석하는 일이 어려워 보일 수도 있다. 한 가지 가능한 전략은 문자열에서 첫 번째 토큰(브라우저의 종류를 어느 정도 알 수 있을 만큼)을 잘라내서 사용자 행동에 대한 또 다른 개요를 만드는 것이다.

results = pd.Series([x.split()[0] for x in frame.a.dropna()])

results[:5]
0               Mozilla/5.0
1    GoogleMaps/RochesterNY
2               Mozilla/4.0
3               Mozilla/5.0
4               Mozilla/5.0
dtype: object
results.value_counts()[:8]
Mozilla/5.0                 2594
Mozilla/4.0                  601
GoogleMaps/RochesterNY       121
Opera/9.80                    34
TEST_INTERNET_AGENT           24
GoogleProducer                21
Mozilla/6.0                    5
BlackBerry8520/5.0.0.681       4
dtype: int64

이제 표준시간대 순위표를 윈도우 사용자와 비윈도우 사용자 그룹으로 나눠보자. 문제를 단순화해서 agent 문자열이 ‘Windows’를 포함하면 윈도우 사용자라고 가정하고 agnet 값이 없는 데이터는 다음과 같이 제외한다.

cframe = frame[frame.a.notnull()]

그리고 이제 각 로우가 윈도우인지 아닌지 검사한다.

import numpy as np

cframe['os'] = np.where(cframe['a'].str.contains('Windows'),
                        'Windows', 'Not Windows')

cframe['os'][:5]
C:\Users\Sangjin\AppData\Local\Temp\ipykernel_19180\543647922.py:3: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
  cframe['os'] = np.where(cframe['a'].str.contains('Windows'),
0        Windows
1    Not Windows
2        Windows
3    Not Windows
4        Windows
Name: os, dtype: object

그런 다음 표준시간대와 운영체제를 기준으로 데이터를 그룹으로 묶는다.

by_tz_os = cframe.groupby(['tz', 'os'])

앞에서 살펴본 value_counts() 함수처럼 그룹별 합계는 size() 함수로 계산할 수 있었다. 결과는 unstack() 함수를 이용해서 표로 재배치한다.

agg_counts = by_tz_os.size().unstack().fillna(0)

agg_counts[:10]
os Not Windows Windows
tz
245.0 276.0
Africa/Cairo 0.0 3.0
Africa/Casablanca 0.0 1.0
Africa/Ceuta 0.0 2.0
Africa/Johannesburg 0.0 1.0
Africa/Lusaka 0.0 1.0
America/Anchorage 4.0 1.0
America/Argentina/Buenos_Aires 1.0 0.0
America/Argentina/Cordoba 0.0 1.0
America/Argentina/Mendoza 0.0 1.0

마지막으로 전체 표준시간대의 순위를 모아보자. 먼저 agg_counts를 보자.

# 오름차순으로 정렬
indexer = agg_counts.sum(1).argsort()

indexer[:10]
tz
                                  24
Africa/Cairo                      20
Africa/Casablanca                 21
Africa/Ceuta                      92
Africa/Johannesburg               87
Africa/Lusaka                     53
America/Anchorage                 54
America/Argentina/Buenos_Aires    57
America/Argentina/Cordoba         26
America/Argentina/Mendoza         55
dtype: int64

agg_counts에 take()를 사용해서 로우를 정렬된 순서 그대로 선택하고 마지막 10개 로우(가장 큰 값)만 잘라낸다.

count_subset = agg_counts.take(indexer[-10:])

count_subset
os Not Windows Windows
tz
America/Sao_Paulo 13.0 20.0
Europe/Madrid 16.0 19.0
Pacific/Honolulu 0.0 36.0
Asia/Tokyo 2.0 35.0
Europe/London 43.0 31.0
America/Denver 132.0 59.0
America/Los_Angeles 130.0 252.0
America/Chicago 115.0 285.0
245.0 276.0
America/New_York 339.0 912.0

pandas에는 이와 똑같은 동작을 하는 nlargest()라는 편리한 메서드가 존재한다.

agg_counts.sum(1).nlargest(10)
tz
America/New_York       1251.0
                        521.0
America/Chicago         400.0
America/Los_Angeles     382.0
America/Denver          191.0
Europe/London            74.0
Asia/Tokyo               37.0
Pacific/Honolulu         36.0
Europe/Madrid            35.0
America/Sao_Paulo        33.0
dtype: float64

그런 다음 앞에서 해본 것처럼 plot() 함수에 stacked=True를 넘겨주면 데이터를 중첩막대그래프로 만들 수 있다(그림 14-2).

# 시각화를 위해 데이터 재배치
count_subset = count_subset.stack()

count_subset.name = 'total'

count_subset = count_subset.reset_index()

count_subset[:10]
tz os total
0 America/Sao_Paulo Not Windows 13.0
1 America/Sao_Paulo Windows 20.0
2 Europe/Madrid Not Windows 16.0
3 Europe/Madrid Windows 19.0
4 Pacific/Honolulu Not Windows 0.0
5 Pacific/Honolulu Windows 36.0
6 Asia/Tokyo Not Windows 2.0
7 Asia/Tokyo Windows 35.0
8 Europe/London Not Windows 43.0
9 Europe/London Windows 31.0
fig = plt.figure(figsize=(15, 10))
ax = fig.add_subplot(1, 1, 1)

sns.barplot(x='total', y='tz', hue='os', data=count_subset, ax=ax)

plt.title('그림 14-2 윈도우 사용자와 비윈도우 사용자별 시간대', size=20, loc='left')
plt.xticks(size=15)
plt.xlabel('mean(total)', size=15)
plt.yticks(size=15)
plt.ylabel('tz', size=15)
plt.legend(fontsize=15)
plt.show()

위 그래프로는 작은 그룹에서 윈도우 사용자의 상대 비율을 확인하기 어렵다. 하지만 각 로우에서 총합을 1로 정규화한 뒤 그래프를 그리면 쉽게 확인할 수 있다.

def norm_total(group):
    group['normed_total'] = group.total / group.total.sum()
    return group

results = count_subset.groupby('tz').apply(norm_total)

정규화한 데이터를 그래프로 그려보자(그림 14-3).

sns.barplot(x='normed_total', y='tz', hue='os', data=results)
plt.xlabel('mean(normed_total)')
plt.show()

groupby()transform() 메서드를 이용해서 정규합 계산을 더 효율적으로 할 수도 있다.

g = count_subset.groupby('tz')

results = count_subset.total / g.total.transform('sum')
results
0     0.393939
1     0.606061
2     0.457143
3     0.542857
4     0.000000
5     1.000000
6     0.054054
7     0.945946
8     0.581081
9     0.418919
10    0.691099
11    0.308901
12    0.340314
13    0.659686
14    0.287500
15    0.712500
16    0.470250
17    0.529750
18    0.270983
19    0.729017
Name: total, dtype: float64

MovieLens의 영화 평점 데이터

GroupLens 연구소는 1990년대 말부터 2000년대 초까지 MovieLens 사용자로부터 수집한 방대한 영화 평점 데이터를 제공하고 있다. 이 데이터에는 영화 평점과 영화에 대한 정보(장르, 개봉 년도) 그리고 사용자에 대한 정보(나이, 우편번호, 성별, 직업)가 포함되어 있다. 이런 종류의 데이터는 머신러닝 알고리즘 기반의 추천 시스템을 개발하는 데 주로 활용한다. 머신러닝 기법을 여기서 소개하기는 어렵고 이런 종류의 데이터를 요구 사항에 맞도록 잘 쪼개는 방법만 소개하겠다.

MovieLens 1M(백만 개) 데이터셋은 약 6,000여 명의 사용자로부터 수집한 4,000여 편의 영화에 대한 백만 개의 영화 평점을 담고 있다. 이 데이터셋은 평점, 사용자 정보, 영화 정보의 3가지 테이블로 나뉘어 있는데, zip 파일의 압축을 풀고 각 테이블을 pandas.read_table() 함수를 사용하여 DataFrame 객체로 불러오자.

import pandas as pd

# 출력되는 내용을 줄인다.
pd.options.display.max_rows = 10

unames = ['user_id', 'gender', 'age', 'occupation', 'zip']
users = pd.read_table('datasets/movielens/users.dat', sep='::',
                      header=None, names=unames, engine='python')

rnames = ['user_id', 'movie_id', 'rating', 'timestamp']
ratings = pd.read_table('datasets/movielens/ratings.dat', sep='::',
                        header=None, names=rnames, engine='python')

mnaems = ['movie_id', 'title', 'genres']
movies = pd.read_table('datasets/movielens/movies.dat', sep='::',
                       header=None, names=mnaems, engine='python')

DataFrame 객체에 데이터가 제대로 들어갔는지 확인하기 위해 파이썬의 리스트 분할 문법을 사용해서 첫 5개 로우를 출력해보자.

users[:5]
user_id gender age occupation zip
0 1 F 1 10 48067
1 2 M 56 16 70072
2 3 M 25 15 55117
3 4 M 45 7 02460
4 5 M 25 20 55455
ratings[:5]
user_id movie_id rating timestamp
0 1 1193 5 978300760
1 1 661 3 978302109
2 1 914 3 978301968
3 1 3408 4 978300275
4 1 2355 5 978824291
movies[:5]
movie_id title genres
0 1 Toy Story (1995) Animation|Children's|Comedy
1 2 Jumanji (1995) Adventure|Children's|Fantasy
2 3 Grumpier Old Men (1995) Comedy|Romance
3 4 Waiting to Exhale (1995) Comedy|Drama
4 5 Father of the Bride Part II (1995) Comedy

나이와 직업은 실제 값이 아니라 그룹을 가리키는 코드 번호이며 데이터셋에 있는 README 파일에 해당 코드에 대한 설명이 들어 있다. 세 종류의 테이블에 걸쳐 있는 데이터를 분석하는 일은 단순한 작업이 아니다. 나이와 성별에 따른 어떤 영화의 평균 평점을 계산한다고 해보자. 아래에서 확인할 수 있지만, 모든 데이터를 하나의 테이블로 병합하여 계산하면 무척 쉽게 처리할 수 있다. pandas의 merge() 함수를 이용해서 ratings 테이블과 users 테이블을 병합하고 그 결과를 다시 movies 테이블과 병합한다. pandas는 병합하려는 두 테이블에서 중복되는 칼럼의 이름을 키로 사용한다.

data = pd.merge(pd.merge(ratings, users), movies)

data
user_id movie_id rating timestamp gender age occupation zip title genres
0 1 1193 5 978300760 F 1 10 48067 One Flew Over the Cuckoo's Nest (1975) Drama
1 2 1193 5 978298413 M 56 16 70072 One Flew Over the Cuckoo's Nest (1975) Drama
2 12 1193 4 978220179 M 25 12 32793 One Flew Over the Cuckoo's Nest (1975) Drama
3 15 1193 4 978199279 M 25 7 22903 One Flew Over the Cuckoo's Nest (1975) Drama
4 17 1193 5 978158471 M 50 1 95350 One Flew Over the Cuckoo's Nest (1975) Drama
... ... ... ... ... ... ... ... ... ... ...
1000204 5949 2198 5 958846401 M 18 17 47901 Modulations (1998) Documentary
1000205 5675 2703 3 976029116 M 35 14 30030 Broken Vessels (1998) Drama
1000206 5780 2845 1 958153068 M 18 17 92886 White Boys (1999) Drama
1000207 5851 3607 5 957756608 F 18 20 55410 One Little Indian (1973) Comedy|Drama|Western
1000208 5938 2909 4 957273353 M 25 1 35401 Five Wives, Three Secretaries and Me (1998) Documentary

1000209 rows × 10 columns

data.iloc[0]
user_id                                            1
movie_id                                        1193
rating                                             5
timestamp                                  978300760
gender                                             F
age                                                1
occupation                                        10
zip                                            48067
title         One Flew Over the Cuckoo's Nest (1975)
genres                                         Drama
Name: 0, dtype: object

성별에 따른 각 영화의 평균 평점을 구하려면 pivot_table() 메서드를 사용하면 된다.

mean_ratings = data.pivot_table('rating', index='title', 
                                columns='gender', aggfunc='mean')

mean_ratings[:5]
gender F M
title
$1,000,000 Duck (1971) 3.375000 2.761905
'Night Mother (1986) 3.388889 3.352941
'Til There Was You (1997) 2.675676 2.733333
'burbs, The (1989) 2.793478 2.962085
...And Justice for All (1979) 3.828571 3.689024

이렇게 하면 매 로우마다 성별에 따른 평균 영화 평점 정보를 담고 있는 DataFrame 객체가 만들어진다. 먼저 250건 이상의 평점 정보가 있는 영화만 추려보자. 데이터를 영화 제목으로 그룹화하고 size() 함수를 사용해서 제목별 평점 정보 건수를 Series 객체로 얻어낸다.

ratings_by_title = data.groupby('title').size()

ratings_by_title[:10]
title
$1,000,000 Duck (1971)                37
'Night Mother (1986)                  70
'Til There Was You (1997)             52
'burbs, The (1989)                   303
...And Justice for All (1979)        199
1-900 (1994)                           2
10 Things I Hate About You (1999)    700
101 Dalmatians (1961)                565
101 Dalmatians (1996)                364
12 Angry Men (1957)                  616
dtype: int64
active_titles = ratings_by_title.index[ratings_by_title >= 250]

active_titles
Index([''burbs, The (1989)', '10 Things I Hate About You (1999)',
       '101 Dalmatians (1961)', '101 Dalmatians (1996)', '12 Angry Men (1957)',
       '13th Warrior, The (1999)', '2 Days in the Valley (1996)',
       '20,000 Leagues Under the Sea (1954)', '2001: A Space Odyssey (1968)',
       '2010 (1984)',
       ...
       'X-Men (2000)', 'Year of Living Dangerously (1982)',
       'Yellow Submarine (1968)', 'You've Got Mail (1998)',
       'Young Frankenstein (1974)', 'Young Guns (1988)',
       'Young Guns II (1990)', 'Young Sherlock Holmes (1985)',
       'Zero Effect (1998)', 'eXistenZ (1999)'],
      dtype='object', name='title', length=1216)

250건 이상의 평점 정보가 있는 영화에 대한 색인은 mean_ratings에서 항목을 선택하기 위해 사용할 수 있다.

# 영화 색인으로 로우 선택
mean_ratings = mean_ratings.loc[active_titles]

mean_ratings
gender F M
title
'burbs, The (1989) 2.793478 2.962085
10 Things I Hate About You (1999) 3.646552 3.311966
101 Dalmatians (1961) 3.791444 3.500000
101 Dalmatians (1996) 3.240000 2.911215
12 Angry Men (1957) 4.184397 4.328421
... ... ...
Young Guns (1988) 3.371795 3.425620
Young Guns II (1990) 2.934783 2.904025
Young Sherlock Holmes (1985) 3.514706 3.363344
Zero Effect (1998) 3.864407 3.723140
eXistenZ (1999) 3.098592 3.289086

1216 rows × 2 columns

여성에게 높은 평점을 받은 영화 목록을 확인하기 위해 F 칼럼을 내림차순으로 정렬한다.

top_female_ratings = mean_ratings.sort_values(by='F', ascending=False)

top_female_ratings[:10]
gender F M
title
Close Shave, A (1995) 4.644444 4.473795
Wrong Trousers, The (1993) 4.588235 4.478261
Sunset Blvd. (a.k.a. Sunset Boulevard) (1950) 4.572650 4.464589
Wallace & Gromit: The Best of Aardman Animation (1996) 4.563107 4.385075
Schindler's List (1993) 4.562602 4.491415
Shawshank Redemption, The (1994) 4.539075 4.560625
Grand Day Out, A (1992) 4.537879 4.293255
To Kill a Mockingbird (1962) 4.536667 4.372611
Creature Comforts (1990) 4.513889 4.272277
Usual Suspects, The (1995) 4.513317 4.518248

평점 차이 구하기

이번엔 남녀 간의 호불호가 갈리는 영화를 찾아보자. mean_ratings에 평균 평점의 차이를 담을 수 있는 칼럼을 하나 추가하고, 그 칼럼을 기준으로 정렬하자.

mean_ratings['diff'] = mean_ratings['M'] - mean_ratings['F']

이제 ‘diff’로 정렬하면 여성들이 더 선호하는 영화 순으로 나타난다.

sorted_by_diff = mean_ratings.sort_values(by='diff')

sorted_by_diff[:10]
gender F M diff
title
Dirty Dancing (1987) 3.790378 2.959596 -0.830782
Jumpin' Jack Flash (1986) 3.254717 2.578358 -0.676359
Grease (1978) 3.975265 3.367041 -0.608224
Little Women (1994) 3.870588 3.321739 -0.548849
Steel Magnolias (1989) 3.901734 3.365957 -0.535777
Anastasia (1997) 3.800000 3.281609 -0.518391
Rocky Horror Picture Show, The (1975) 3.673016 3.160131 -0.512885
Color Purple, The (1985) 4.158192 3.659341 -0.498851
Age of Innocence, The (1993) 3.827068 3.339506 -0.487561
Free Willy (1993) 2.921348 2.438776 -0.482573

역순으로 정렬한 다음 상위 10개 로우를 잘라내면 남성의 선호도 순으로도 확인할 수 있다.

# 뒤에서 10개의 로우 선택
sorted_by_diff[::-1][:10]
gender F M diff
title
Good, The Bad and The Ugly, The (1966) 3.494949 4.221300 0.726351
Kentucky Fried Movie, The (1977) 2.878788 3.555147 0.676359
Dumb & Dumber (1994) 2.697987 3.336595 0.638608
Longest Day, The (1962) 3.411765 4.031447 0.619682
Cable Guy, The (1996) 2.250000 2.863787 0.613787
Evil Dead II (Dead By Dawn) (1987) 3.297297 3.909283 0.611985
Hidden, The (1987) 3.137931 3.745098 0.607167
Rocky III (1982) 2.361702 2.943503 0.581801
Caddyshack (1980) 3.396135 3.969737 0.573602
For a Few Dollars More (1965) 3.409091 3.953795 0.544704

성별에 관계없이 영화에 대한 호불호가 극명하게 나뉘는 영화를 찾아보자. 호불호는 평점의 분산이나 표준편차로 측정할 수 있다.

# 영화별 평점 표준편차
rating_std_by_title = data.groupby('title')['rating'].std()

# active_titles만 선택
rating_std_by_title = rating_std_by_title.loc[active_titles]

# 평점 내림차순으로 Series 정렬
rating_std_by_title.sort_values(ascending=False)[:10]
title
Dumb & Dumber (1994)                     1.321333
Blair Witch Project, The (1999)          1.316368
Natural Born Killers (1994)              1.307198
Tank Girl (1995)                         1.277695
Rocky Horror Picture Show, The (1975)    1.260177
Eyes Wide Shut (1999)                    1.259624
Evita (1996)                             1.253631
Billy Madison (1995)                     1.249970
Fear and Loathing in Las Vegas (1998)    1.246408
Bicentennial Man (1999)                  1.245533
Name: rating, dtype: float64

데이터 파일을 열어본 독자라면 영화 장르가 버티컬 바(|)로 구분되어 제공되고 있다는 사실을 알 수 있다. 만일 영화 장르에 기반한 분석을 하려면 영화 장르 정보를 좀 더 사용하기 편한 형식으로 변형할 필요가 있다.

신생아 이름

미국사회보장국(SSA)에서는 1880년부터 현재까지 가장 빈도가 높은 신생아 이름에 대한 정보를 제공하고 있다. 여러 가지 유명한 R 패키지 개발자인 해들리 위컴은 R에서 데이터를 다루는 방법을 설명할 때 종종 이 데이터셋을 활용한다.

이 데이터셋을 불러 오려면 데이터를 다듬는 과정이 필요한데 일단 정리하고 나면 다음과 같은 DataFrame을 얻을 수 있다.

names.head(10)
name sex births year
0 Mary F 7065 1880
1 Anna F 2604 1880
2 Emma F 2003 1880
3 Elizabeth F 1939 1880
4 Minnie F 1746 1880
5 Margaret F 1578 1880
6 Ida F 1472 1880
7 Alice F 1414 1880
8 Bertha F 1320 1880
9 Sarah F 1288 1880

이 데이터를 이용해서 여러 가지 분석을 할 수 있다.

  • 시대별로 특정 이름이 차지하는 비율을 구해 얼마나 흔한 이름인지 알아보기

  • 이름의 상대 순위 알아보기

  • 각 연도별로 가장 인기 있는 이름, 가장 많이 증가하거나 감소한 이름 알아보기

  • 모음, 자음, 길이, 전체 다양성, 철자 변화, 첫 글자와 마지막 글자 등 이름 유행 분석하기

  • 성서에 등장하는 이름, 유명인, 인구통계학적 변화 등 외부 자료를 통한 유행 분석

지금까지 살펴본 도구를 이용하면 이 정도의 분석은 아주 쉽게 해낼 수 있다. 여러분이 따라 할 수 있게 차근차근 설명하겠다.

저자가 이 책을 쓰는 시기 미국사회보장국은 매해 전체 출생에 대해 성별과 이름 정보를 제공하고 있다. 가공되지 않은 데이터 파일은 사회보장국 웹페이지1에서 구할 수 있다.

페이지 주소는 바뀔 수 있으니 인터넷 검색을 통해 찾아보자. ‘국가 정보’인 names.zip 파일을 받은 후 압축을 풀면 yob1880.txt 같은 이름의 파일들이 담겨 있는 디렉터리가 생성된다. 유닉스의 head 명령어를 사용하여 이 파일들 중 하나에서 처음 10줄을 먼저 살펴보자(윈도우에서는 more 명령어를 사용하거나 텍스트 데이터에서 열어본다).

%load datasets/babynames/yob1880.txt

Mary,F,7065
Anna,F,2604
Emma,F,2003
Elizabeth,F,1939
Minnie,F,1746

위 매직커맨드를 실행시켜 확인해보면 파일이 편리하게도 쉼표로 구분되어 있으니 pandas.read_csv() 메서드를 사용해서 DataFrame 객체로 불러오자.

import pandas as pd

names1880 = pd.read_csv('datasets/babynames/yob1880.txt',
                        names=['name', 'sex', 'births'])

names1880
name sex births
0 Mary F 7065
1 Anna F 2604
2 Emma F 2003
3 Elizabeth F 1939
4 Minnie F 1746
... ... ... ...
1995 Woodie M 5
1996 Worthy M 5
1997 Wright M 5
1998 York M 5
1999 Zachariah M 5

2000 rows × 3 columns

이 데이터는 각 연도별로 최소 5명 이상 중복되는 이름만 포함하고 있다. 따라서 편의상 성별별 출생수를 모두 합한 값을 해당 연도의 전체 출생수라고 가정하자.

names1880.groupby('sex').births.sum()
sex
F     90993
M    110493
Name: births, dtype: int64

자료가 연도별 파일로 나뉘어져 있으니 먼저 모든 데이터를 DataFrame 하나로 모은 다음 year 항목을 추가한다. pandas.concat()을 이용하면 이 작업을 쉽게 할 수 있다.

years = range(1880, 2011)

pieces = []
columns = ['name', 'sex', 'births']

for year in years:
    path = 'datasets/babynames/yob%d.txt' % year
    frame = pd.read_csv(path, names=columns)
    
    frame['year'] = year
    pieces.append(frame)
    
# 모두 하나의 DataFrame으로 합치기
names = pd.concat(pieces, ignore_index=True)

여기서 두 가지 언급해야 할 내용이 있다. 첫째, concat() 메서드는 DataFrame 객체를 합쳐준다.

둘째, read_csv()로 읽어온 원래 로우 순서는 몰라도 되니 concat() 메서드에 ignore_index=True를 인자로 전달해야 한다. 이렇게 해서 전체 이름 데이터를 담고 있는 거대한 DataFrame 객체를 만들었다.

names
name sex births year
0 Mary F 7065 1880
1 Anna F 2604 1880
2 Emma F 2003 1880
3 Elizabeth F 1939 1880
4 Minnie F 1746 1880
... ... ... ... ...
1690779 Zymaire M 5 2010
1690780 Zyonne M 5 2010
1690781 Zyquarius M 5 2010
1690782 Zyran M 5 2010
1690783 Zzyzx M 5 2010

1690784 rows × 4 columns

이제 이 데이터에 groupby()pivot_table()을 이용해서 연도나 성별에 따른 데이터를 수집할 수 있다(그림 14-4).

total_births = names.pivot_table('births', index='year',
                                 columns='sex', aggfunc=sum)

total_births.tail()
sex F M
year
2006 1896468 2050234
2007 1916888 2069242
2008 1883645 2032310
2009 1827643 1973359
2010 1759010 1898382

그림 14-4 연도와 성별별 출산수

total_births.plot(title='Total births by sex and year')
print('그림 14-4 연도와 성별별 출산수')
그림 14-4 연도와 성별별 출산수

다음은 prop 칼럼을 추가해서 각 이름이 전체 출생수에서 차지하는 비율을 계산하자. prop 값이 0.02라면 100명의 아기 중 2명의 이름이 같다는 뜻이다. 데이터를 연도와 성별로 그룹화하고 각 그룹에 새 칼럼을 추가하자.

def add_prop(group):
    group['prop'] = group.births / group.births.sum()
    return group

names = names.groupby(['year', 'sex']).apply(add_prop)

이제 names에 새로운 prop 칼럼이 추가되었다.

names
name sex births year prop
0 Mary F 7065 1880 0.077643
1 Anna F 2604 1880 0.028618
2 Emma F 2003 1880 0.022013
3 Elizabeth F 1939 1880 0.021309
4 Minnie F 1746 1880 0.019188
... ... ... ... ... ...
1690779 Zymaire M 5 2010 0.000003
1690780 Zyonne M 5 2010 0.000003
1690781 Zyquarius M 5 2010 0.000003
1690782 Zyran M 5 2010 0.000003
1690783 Zzyzx M 5 2010 0.000003

1690784 rows × 5 columns

그룹 관련 연산을 수행할 때는 모든 그룹에서 prop 칼럼의 합이 1이 맞는지 확인하는 새너티 테스트를 하는 게 좋다.

names.groupby(['year', 'sex']).prop.sum()
year  sex
1880  F      1.0
      M      1.0
1881  F      1.0
      M      1.0
1882  F      1.0
            ... 
2008  M      1.0
2009  F      1.0
      M      1.0
2010  F      1.0
      M      1.0
Name: prop, Length: 262, dtype: float64

이제 모든 준비가 끝났고, 분석에 사용할 각 연도별/성별에 따른 선호하는 이름 1,000개를 추출하자. 이것 역시 그룹 연산이다.

def get_top1000(group):
    return group.sort_values(by='births', ascending=False)[:1000]

grouped = names.groupby(['year', 'sex'])
top1000 = grouped.apply(get_top1000)
# 그룹 색인은 필요없으므로 삭제
top1000.reset_index(inplace=True, drop=True)

함수를 정의하지 않고 직접 추출하고 싶다면 다음처럼 할 수도 있다.

pieces = []
for year, group in names.groupby(['year', 'sex']):
    pieces.append(group.sort_values(by='births', ascending=False)[:1000])
top1000 = pd.concat(pieces, ignore_index=True)

이제 데이터셋의 크기가 조금 줄었다.

top1000
name sex births year prop
0 Mary F 7065 1880 0.077643
1 Anna F 2604 1880 0.028618
2 Emma F 2003 1880 0.022013
3 Elizabeth F 1939 1880 0.021309
4 Minnie F 1746 1880 0.019188
... ... ... ... ... ...
261872 Camilo M 194 2010 0.000102
261873 Destin M 194 2010 0.000102
261874 Jaquan M 194 2010 0.000102
261875 Jaydan M 194 2010 0.000102
261876 Maxton M 193 2010 0.000102

261877 rows × 5 columns

이렇게 추출한 상위 1,000개의 이름 데이터는 이어지는 분석에서 사용하도록 하자.

이름 유행 분석

전체 데이터셋과 상위 1,000개의 이름 데이터로 흥미로운 이름 유행을 분석해보자. 먼저 상위 1,000개의 데이터를 남자아이와 여자아이로 분리한다.

boys = top1000[top1000.sex == 'M']

girs = top1000[top1000.sex == 'F']

연도별로 John이나 Mary라는 이름의 추이를 간단하게 그래프로 그릴 수 있는데, 그전에 데이터를 살짝 변경할 필요가 있다. 연도와 이름에 대한 전체 출생수를 피벗테이블로 만들자.

total_births = top1000.pivot_table('births', index='year',
                                   columns='name',
                                   aggfunc=sum)

total_births
name Aaden Aaliyah Aarav Aaron Aarush Ab Abagail Abb Abbey Abbie ... Zoa Zoe Zoey Zoie Zola Zollie Zona Zora Zula Zuri
year
1880 NaN NaN NaN 102.0 NaN NaN NaN NaN NaN 71.0 ... 8.0 23.0 NaN NaN 7.0 NaN 8.0 28.0 27.0 NaN
1881 NaN NaN NaN 94.0 NaN NaN NaN NaN NaN 81.0 ... NaN 22.0 NaN NaN 10.0 NaN 9.0 21.0 27.0 NaN
1882 NaN NaN NaN 85.0 NaN NaN NaN NaN NaN 80.0 ... 8.0 25.0 NaN NaN 9.0 NaN 17.0 32.0 21.0 NaN
1883 NaN NaN NaN 105.0 NaN NaN NaN NaN NaN 79.0 ... NaN 23.0 NaN NaN 10.0 NaN 11.0 35.0 25.0 NaN
1884 NaN NaN NaN 97.0 NaN NaN NaN NaN NaN 98.0 ... 13.0 31.0 NaN NaN 14.0 6.0 8.0 58.0 27.0 NaN
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
2006 NaN 3737.0 NaN 8279.0 NaN NaN 297.0 NaN 404.0 440.0 ... NaN 5145.0 2839.0 530.0 NaN NaN NaN NaN NaN NaN
2007 NaN 3941.0 NaN 8914.0 NaN NaN 313.0 NaN 349.0 468.0 ... NaN 4925.0 3028.0 526.0 NaN NaN NaN NaN NaN NaN
2008 955.0 4028.0 219.0 8511.0 NaN NaN 317.0 NaN 344.0 400.0 ... NaN 4764.0 3438.0 492.0 NaN NaN NaN NaN NaN NaN
2009 1265.0 4352.0 270.0 7936.0 NaN NaN 296.0 NaN 307.0 369.0 ... NaN 5120.0 3981.0 496.0 NaN NaN NaN NaN NaN NaN
2010 448.0 4628.0 438.0 7374.0 226.0 NaN 277.0 NaN 295.0 324.0 ... NaN 6200.0 5164.0 504.0 NaN NaN NaN NaN NaN 258.0

131 rows × 6868 columns

DataFrame의 plot() 메서드를 사용해서 몇몇 이름의 추이를 그래프로 그려보자(그림 14-5).

total_births.info()
<class 'pandas.core.frame.DataFrame'>
Int64Index: 131 entries, 1880 to 2010
Columns: 6868 entries, Aaden to Zuri
dtypes: float64(6868)
memory usage: 6.9 MB
subset = total_births[['John', 'Harry', 'Mary', 'Marilyn']]

subset.plot(subplots=True, figsize=(12, 10), grid=False,
            title="Number of births per year")
print('그림 14-5 연도별 남아 및 여아 이름 추이')
그림 14-5 연도별 남아 및 여아 이름 추이

그래프를 보면 예로 든 이름들이 최근 미국에서 인기가 없다는 걸 알 수 있다. 하지만 단순히 이렇게 결론을 내기에는 조금 복잡한데 이 부분은 다음 절에서 살펴보겠다.

다양한 이름을 사용하는 경향 측정하기

위에서 확인한 그래프의 감소 추세는 부모가 아이의 이름을 지을 때 흔한 이름은 기피한다고 해석할 수 있다. 이 가설은 데이터에서 살펴볼 수 있으며 확인이 가능하다. 좀 더 자세히 알아보기 위해 인기 있는 이름 1,000개가 전체 출생수에서 차지하는 비율을 연도별/성별 그래프로 그려보자(그림 14-6).

table = top1000.pivot_table('prop', index='year',
                            columns='sex', aggfunc=sum)

table
sex F M
year
1880 1.000000 0.997375
1881 1.000000 1.000000
1882 0.998702 0.995646
1883 0.997596 0.998566
1884 0.993156 0.994539
... ... ...
2006 0.753153 0.860368
2007 0.745959 0.855159
2008 0.740933 0.850003
2009 0.737290 0.845256
2010 0.736780 0.843156

131 rows × 2 columns

table.plot(title='Sum of table1000.prop by year and sex',
            yticks=np.linspace(0, 1.2, 13), xticks=range(1880, 2020, 10))
print('그림 14-6 인기 있는 이름 1000개의 연도별/성별 비율')
그림 14-6 인기 있는 이름 1000개의 연도별/성별 비율

그림에서 확인할 수 있듯이 실제로 이름의 다양성이 증가하고 있음을 보여준다(상위 1,000개의 이름에서 비율의 총합이 감소하고 있다). 또한 인기 있는 이름순으로 정렬했을 때 전체 출생수의 50%를 차지하기까지 등장하는 이름수도 흥미롭다.

df = boys[boys.year == 2010]

df
name sex births year prop
260877 Jacob M 21875 2010 0.011523
260878 Ethan M 17866 2010 0.009411
260879 Michael M 17133 2010 0.009025
260880 Jayden M 17030 2010 0.008971
260881 William M 16870 2010 0.008887
... ... ... ... ... ...
261872 Camilo M 194 2010 0.000102
261873 Destin M 194 2010 0.000102
261874 Jaquan M 194 2010 0.000102
261875 Jaydan M 194 2010 0.000102
261876 Maxton M 193 2010 0.000102

1000 rows × 5 columns

prop을 내림차순으로 정렬하고 나서 전체의 50%가 되기까지 얼마나 많은 이름이 등장하는지 알아보자. for 문을 사용해서 구현할 수도 있지만, 벡터하된 NumPy를 사용하는 편이 조금 더 편하다. prop의 누계를 cumsum에 저장하고 searchsorted() 메서드를 호출해서 정렬된 상태에서 누계가 0.5가 되는 위치를 구한다.

prop_cumsum = df.sort_values(by='prop', ascending=False).prop.cumsum()

prop_cumsum[:10]
260877    0.011523
260878    0.020934
260879    0.029959
260880    0.038930
260881    0.047817
260882    0.056579
260883    0.065155
260884    0.073414
260885    0.081528
260886    0.089621
Name: prop, dtype: float64
prop_cumsum.values.searchsorted(0.5)
116

배열의 색인은 0부터 시작하기 때문에 결과에 1을 더해주면 117이 나온다. 1900년에는 이보다 더 낮았다.

df = boys[boys.year == 1900]

in1900 = df.sort_values(by='prop', ascending=False).prop.cumsum()

in1900.values.searchsorted(0.5) + 1
25

이제 이 연산을 각 연도별/성별 조합에 적용할 수 있다. 연도와 성을 groupby()로 묶고 각 그룹에 apply()를 사용해서 이 연산을 적용하면 된다.

def get_quantile_count(group, q=0.5):
    group = group.sort_values(by='prop', ascending=False)
    return group.prop.cumsum().values.searchsorted(q) + 1

diversity = top1000.groupby(['year', 'sex']).apply(get_quantile_count)
diversity = diversity.unstack('sex')

연산 결과인 diversity DataFrame은 이제 각 성별에 따라 연도별로 색인된 두 개의 시계열 데이터를 담고 있다.

diversity.head()
sex F M
year
1880 38 14
1881 38 14
1882 38 15
1883 39 15
1884 39 16
diversity.plot(title="Number of popular names in top 50%")
print('그림 14-7 연도별 이름 다양성 지수')
그림 14-7 연도별 이름 다양성 지수

보다시피 여자아이의 이름은 항상 남자아이 이름보다 더 다양하며, 시간이 흐를수록 더욱 다양해지고 있다. 대체되는 철자의 증가 같은 다양성을 높이는 요건에 대한 자세한 분석은 여러분의 몫으로 남겨두겠다.

‘마지막 글자’의 변화

2007년 아이 이름을 연구하는 로라 와튼버그는 지난 100년 동안 남자아이 이름의 마지막 글자의 분포에 중요한 변화가 있었다고 자신의 웹사이트[^2]에 게재했다. 지금부터 전체 자료에서 연도, 성별, 이름의 마지막 글자를 수집해서 이를 확인해보자.

# name 칼럼에서 마지막 글자를 추출한다.
get_last_letter = lambda x: x[-1]
last_letters = names.name.map(get_last_letter)
last_letters.name = 'last_letter'

table = names.pivot_table('births', index=last_letters,
                          columns=['sex', 'year'], aggfunc=sum)

table
sex F ... M
year 1880 1881 1882 1883 1884 1885 1886 1887 1888 1889 ... 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010
last_letter
a 31446.0 31581.0 36536.0 38330.0 43680.0 45408.0 49100.0 48942.0 59442.0 58631.0 ... 39124.0 38815.0 37825.0 38650.0 36838.0 36156.0 34654.0 32901.0 31430.0 28438.0
b NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN ... 50950.0 49284.0 48065.0 45914.0 43144.0 42600.0 42123.0 39945.0 38862.0 38859.0
c NaN NaN 5.0 5.0 NaN NaN NaN NaN NaN NaN ... 27113.0 27238.0 27697.0 26778.0 26078.0 26635.0 26864.0 25318.0 24048.0 23125.0
d 609.0 607.0 734.0 810.0 916.0 862.0 1007.0 1027.0 1298.0 1374.0 ... 60838.0 55829.0 53391.0 51754.0 50670.0 51410.0 50595.0 47910.0 46172.0 44398.0
e 33378.0 34080.0 40399.0 41914.0 48089.0 49616.0 53884.0 54353.0 66750.0 66663.0 ... 145395.0 144651.0 144769.0 142098.0 141123.0 142999.0 143698.0 140966.0 135496.0 129012.0
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
v NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN ... 1209.0 1332.0 1652.0 1823.0 1794.0 2010.0 2295.0 2418.0 2589.0 2723.0
w NaN 5.0 NaN NaN NaN NaN 5.0 NaN NaN NaN ... 52265.0 50103.0 49079.0 47556.0 45464.0 43217.0 40251.0 36937.0 33181.0 30656.0
x NaN NaN NaN 7.0 NaN NaN NaN NaN NaN NaN ... 10691.0 11009.0 11718.0 12399.0 13025.0 13992.0 14306.0 14834.0 16640.0 16352.0
y 10469.0 10404.0 12145.0 12063.0 13917.0 13927.0 14936.0 14980.0 17931.0 17601.0 ... 139109.0 134557.0 130569.0 128367.0 125190.0 123707.0 123397.0 122633.0 112922.0 110425.0
z 106.0 95.0 106.0 141.0 148.0 150.0 202.0 188.0 238.0 277.0 ... 2840.0 2737.0 2722.0 2710.0 2903.0 3086.0 3301.0 3473.0 3633.0 3476.0

26 rows × 262 columns

이제 전체 기간 중 세 지점을 골라 이름의 마지막 글자 몇 개를 출력해보자. 1910, 1960, 2010을 골랐다.

subtable = table.reindex(columns=[1910, 1960, 2010], level='year')

subtable.head()
sex F M
year 1910 1960 2010 1910 1960 2010
last_letter
a 108376.0 691247.0 670605.0 977.0 5204.0 28438.0
b NaN 694.0 450.0 411.0 3912.0 38859.0
c 5.0 49.0 946.0 482.0 15476.0 23125.0
d 6750.0 3729.0 2607.0 22111.0 262112.0 44398.0
e 133569.0 435013.0 313833.0 28655.0 178823.0 129012.0

그다음에는 전체 출생수에서 성별별로 각각의 마지막 글자가 차지하는 비율을 계산하기 위해 전체 출생수로 정규화하자.

subtable.sum()
sex  year
F    1910     396416.0
     1960    2022062.0
     2010    1759010.0
M    1910     194198.0
     1960    2132588.0
     2010    1898382.0
dtype: float64
letter_prop = subtable / subtable.sum()

letter_prop
sex F M
year 1910 1960 2010 1910 1960 2010
last_letter
a 0.273390 0.341853 0.381240 0.005031 0.002440 0.014980
b NaN 0.000343 0.000256 0.002116 0.001834 0.020470
c 0.000013 0.000024 0.000538 0.002482 0.007257 0.012181
d 0.017028 0.001844 0.001482 0.113858 0.122908 0.023387
e 0.336941 0.215133 0.178415 0.147556 0.083853 0.067959
... ... ... ... ... ... ...
v NaN 0.000060 0.000117 0.000113 0.000037 0.001434
w 0.000020 0.000031 0.001182 0.006329 0.007711 0.016148
x 0.000015 0.000037 0.000727 0.003965 0.001851 0.008614
y 0.110972 0.152569 0.116828 0.077349 0.160987 0.058168
z 0.002439 0.000659 0.000704 0.000170 0.000184 0.001831

26 rows × 6 columns

이렇게 구한 이름의 마지막 글자 비율로 성별과 출생 연도에 대한 막대그래프를 그려보자(그림 14-8).

import matplotlib.pyplot as plt

fig, axes = plt.subplots(2, 1, figsize=(10, 8), constrained_layout=True)
letter_prop['M'].plot(kind='bar', rot=0, ax=axes[0], title='Male')
letter_prop['F'].plot(kind='bar', rot=0, ax=axes[1], title='Female',
                      legend=False)
print('그림 14-8 이름의 끝 글자 비율')
그림 14-8 이름의 끝 글자 비율

그래프에서 확인할 수 있듯이 ‘n’으로 끝나는 남자아이 이름의 빈도가 1960년도 이후에 급격하게 증가했다. 이제 세 지점이 아닌 전체 자료에 대해 출생연도와 성별, 남자아이 이름에서 몇 가지 글자로 정규화하고 시계열 데이터로 변환하자.

letter_prop = table / table.sum()

dny_ts = letter_prop.loc[['d', 'n', 'y'], 'M'].T

dny_ts.head()
last_letter d n y
year
1880 0.083055 0.153213 0.075760
1881 0.083247 0.153214 0.077451
1882 0.085340 0.149560 0.077537
1883 0.084066 0.151646 0.079144
1884 0.086120 0.149915 0.080405

이 시계열 데이터를 plot() 메서드를 이용해서 연도별 그래프로 만들어보자(그림 14-9).

dny_ts.plot()
print('그림 14-9 d/n/y로 끝나는 이름을 가진 남아의 연도별 출생 비율')
그림 14-9 d/n/y로 끝나는 이름을 가진 남아의 연도별 출생 비율

남자 이름과 여자 이름이 바뀐 경우

또 다른 재미있는 경향은 예전에는 남자 이름으로 선호되다가 현재는 여자 이름으로 선호되는 경우다. Lesley 또는 Leslie 라는 이름이 그렇다. top1000 데이터를 이용해서 ‘lesl’로 시작하는 이름을 포함하는 목록을 만들어보자.

all_names = pd.Series(top1000.name.unique())

lesley_like = all_names[all_names.str.lower().str.contains('lesl')]

lesley_like
632     Leslie
2294    Lesley
4262    Leslee
4728     Lesli
6103     Lesly
dtype: object

이제 이 이름들만 걸러내서 이름별로 출생수를 구하고 상대도수를 확인해보자.

filtered = top1000[top1000.name.isin(lesley_like)]

filtered.groupby('name').births.sum()
name
Leslee      1082
Lesley     35022
Lesli        929
Leslie    370429
Lesly      10067
Name: births, dtype: int64

그리고 성별과 연도별로 모은 다음 출생연도로 정규화한다.

table = filtered.pivot_table('births', index='year',
                             columns='sex', aggfunc='sum')

table = table.div(table.sum(1), axis=0)

table.tail()
sex F M
year
2006 1.0 NaN
2007 1.0 NaN
2008 1.0 NaN
2009 1.0 NaN
2010 1.0 NaN

마지막으로 시대별로 성별에 따른 명세를 그래프로 그려보자(그림 14-10).

table.plot(style={'M': 'k-', 'F': 'k--'})
print('그림 14-10 Lesley와 비슷한 이름의 비율')
그림 14-10 Lesley와 비슷한 이름의 비율

미국농무부 영양소 정보

미국농무부USDA는 음식의 영양소 정보 데이터베이스를 제공하고 있다. 개발자인 애슐리 윌리엄스는 이 데이터베이스를 다음과 같은 JSON 형식으로 제공하고 있다.

db[0]
{'id': 1008,
 'description': 'Cheese, caraway',
 'tags': [],
 'manufacturer': '',
 'group': 'Dairy and Egg Products',
 'portions': [{'amount': 1, 'unit': 'oz', 'grams': 28.35}],
 'nutrients': [{'value': 25.18,
   'units': 'g',
   'description': 'Protein',
   'group': 'Composition'
   },
   ...
  ]
}

각 음식은 숫자로 된 고유 ID뿐만 아니라 영양소와 제공량 등 두 가지 리스트를 가지고 있다. 이러한 형식의 데이터는 분석하기 편하지 않으므로 데이터의 모양을 좀 더 나은 형태로 바꿔보자.

웹사이트에서 데이터를 내려받은 다음 압축을 풀고 선호하는 JSON 라이브러리를 사용해서 파이썬에서 읽어오자. 나는 파이썬에 내장된 json 모듈을 사용하겠다.

import json

db = json.load(open('datasets/usda_food/database.json'))

len(db)
6636

db에 있는 각 엔트리는 한 가지 음식에 대한 모든 정보를 담고 있는 사전형이다. ‘nutrients’ 필드는 사전의 리스트이며 각 항목은 한 가지 영양소에 대한 정보를 담고 있다.

db[0].keys()
dict_keys(['id', 'description', 'tags', 'manufacturer', 'group', 'portions', 'nutrients'])
db[0]['nutrients'][0]
{'value': 25.18,
 'units': 'g',
 'description': 'Protein',
 'group': 'Composition'}
nutrients = pd.DataFrame(db[0]['nutrients'])

nutrients[:7]
value units description group
0 25.18 g Protein Composition
1 29.20 g Total lipid (fat) Composition
2 3.06 g Carbohydrate, by difference Composition
3 3.28 g Ash Other
4 376.00 kcal Energy Energy
5 39.28 g Water Composition
6 1573.00 kJ Energy Energy

사전의 리스트를 DataFrame으로 바꿀 때 추출할 필드 목록을 지정해줄 수 있다. 우리는 음식의 이름과 그룹, id 그리고 제조사를 추출하겠다.

info_keys = ['description', 'group', 'id', 'manufacturer']

info = pd.DataFrame(db, columns=info_keys)

info[:5]
description group id manufacturer
0 Cheese, caraway Dairy and Egg Products 1008
1 Cheese, cheddar Dairy and Egg Products 1009
2 Cheese, edam Dairy and Egg Products 1018
3 Cheese, feta Dairy and Egg Products 1019
4 Cheese, mozzarella, part skim milk Dairy and Egg Products 1028
info.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 6636 entries, 0 to 6635
Data columns (total 4 columns):
 #   Column        Non-Null Count  Dtype 
---  ------        --------------  ----- 
 0   description   6636 non-null   object
 1   group         6636 non-null   object
 2   id            6636 non-null   int64 
 3   manufacturer  5195 non-null   object
dtypes: int64(1), object(3)
memory usage: 207.5+ KB

value_counts() 메서드를 이용해서 음식 그룹의 분포를 확인할 수 있다.

pd.value_counts(info.group)[:10]
Vegetables and Vegetable Products    812
Beef Products                        618
Baked Products                       496
Breakfast Cereals                    403
Legumes and Legume Products          365
Fast Foods                           365
Lamb, Veal, and Game Products        345
Sweets                               341
Fruits and Fruit Juices              328
Pork Products                        328
Name: group, dtype: int64

모든 영양소 정보를 분석해보자. 좀 더 쉬운 분석을 위해 각 음식의 영양소 정보를 거대한 테이블 하나에 담아보자. 그러려면 사전에 몇 가지 과정을 거쳐야 하는데 먼저 음식의 영양소 리스트를 하나의 DataFrame으로 변환하고, 음식의 id를 위한 칼럼을 하나 추가한다. 또한 이 DataFrame을 리스트에 추가한다. 그리고 마지막으로 이 리스트를 concat() 메서드를 사용해서 하나로 합친다.

문제가 없다면 nutrients는 다음과 같은 모습일 것이다.

pieces = []

for i in range(len(db)):
    iid = db[i]['id']
    df = pd.DataFrame(db[i]['nutrients'], columns=['description', 'group', 'units', 'value'])
    df['id'] = iid
    pieces.append(df)
    
# 모두 하나의 DataFrame으로 합치기
nutrients = pd.concat(pieces, ignore_index=True)

nutrients
description group units value id
0 Protein Composition g 25.180 1008
1 Total lipid (fat) Composition g 29.200 1008
2 Carbohydrate, by difference Composition g 3.060 1008
3 Ash Other g 3.280 1008
4 Energy Energy kcal 376.000 1008
... ... ... ... ... ...
389350 Vitamin B-12, added Vitamins mcg 0.000 43546
389351 Cholesterol Other mg 0.000 43546
389352 Fatty acids, total saturated Other g 0.072 43546
389353 Fatty acids, total monounsaturated Other g 0.028 43546
389354 Fatty acids, total polyunsaturated Other g 0.041 43546

389355 rows × 5 columns

이유야 어쨌든 이 DataFrame에는 붕복된 데이터가 있으므로 제거하도록 한다.

nutrients.duplicated().sum()  # 중복 확인
14179
nutrients = nutrients.drop_duplicates()

‘group’과 ‘description’은 모두 DataFrame 객체이므로 뭐가 뭔지 쉽게 알아볼 수 있도록 이름을 바꿔주자.

col_mapping = {'description': 'food',
               'group': 'fgroup'}

info = info.rename(columns=col_mapping, copy=False)

info.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 6636 entries, 0 to 6635
Data columns (total 4 columns):
 #   Column        Non-Null Count  Dtype 
---  ------        --------------  ----- 
 0   food          6636 non-null   object
 1   fgroup        6636 non-null   object
 2   id            6636 non-null   int64 
 3   manufacturer  5195 non-null   object
dtypes: int64(1), object(3)
memory usage: 207.5+ KB
col_mapping = {'description': 'nutrient',
               'group': 'nutgroup'}

nutrients = nutrients.rename(columns=col_mapping, copy=False)

nutrients
nutrient nutgroup units value id
0 Protein Composition g 25.180 1008
1 Total lipid (fat) Composition g 29.200 1008
2 Carbohydrate, by difference Composition g 3.060 1008
3 Ash Other g 3.280 1008
4 Energy Energy kcal 376.000 1008
... ... ... ... ... ...
389350 Vitamin B-12, added Vitamins mcg 0.000 43546
389351 Cholesterol Other mg 0.000 43546
389352 Fatty acids, total saturated Other g 0.072 43546
389353 Fatty acids, total monounsaturated Other g 0.028 43546
389354 Fatty acids, total polyunsaturated Other g 0.041 43546

375176 rows × 5 columns

여기까지 했으면 info 객체를 nutrients 객체로 병합하자.

ndata = pd.merge(nutrients, info, on='id', how='outer')

ndata.info()
<class 'pandas.core.frame.DataFrame'>
Int64Index: 375176 entries, 0 to 375175
Data columns (total 8 columns):
 #   Column        Non-Null Count   Dtype  
---  ------        --------------   -----  
 0   nutrient      375176 non-null  object 
 1   nutgroup      375176 non-null  object 
 2   units         375176 non-null  object 
 3   value         375176 non-null  float64
 4   id            375176 non-null  int64  
 5   food          375176 non-null  object 
 6   fgroup        375176 non-null  object 
 7   manufacturer  293054 non-null  object 
dtypes: float64(1), int64(1), object(6)
memory usage: 25.8+ MB
ndata.iloc[30000]
nutrient                                       Glycine
nutgroup                                   Amino Acids
units                                                g
value                                             0.04
id                                                6158
food            Soup, tomato bisque, canned, condensed
fgroup                      Soups, Sauces, and Gravies
manufacturer                                          
Name: 30000, dtype: object

이제 음식 그룹과 영양소 종류별 중간값을 그래프로 만들 수 있다(그림 14-11).

result = ndata.groupby(['nutrient', 'fgroup'])['value'].quantile(0.5)

result['Zinc, Zn'].sort_values().plot(kind='barh')

print('그림 14-11 음식 그룹별 아연 함량의 중간값')
그림 14-11 음식 그룹별 아연 함량의 중간값

좀 더 응용하면 각 영양소가 어떤 음식에 가장 많이 들어 있는지 찾아볼 수도 있다.

by_nutrient = ndata.groupby(['nutgroup', 'nutrient'])

get_maximum = lambda x: x.loc[x.value.idxmax()]
get_minimum = lambda x: x.loc[x.value.idxmin()]

max_foods = by_nutrient.apply(get_maximum)[['value', 'food']]

# 음식의 종류를 제한하자.
max_foods.food = max_foods.food.str[:50]

위 코드의 결과 DataFrame을 이 부분에 다 싣기에는 너무 방대하므로 아미노산Amino Acids에 대한 내용만 싣겠다.

max_foods.loc['Amino Acids']['food']
nutrient
Alanine                          Gelatins, dry powder, unsweetened
Arginine                              Seeds, sesame flour, low-fat
Aspartic acid                                  Soy protein isolate
Cystine               Seeds, cottonseed flour, low fat (glandless)
Glutamic acid                                  Soy protein isolate
                                       ...                        
Serine           Soy protein isolate, PROTEIN TECHNOLOGIES INTE...
Threonine        Soy protein isolate, PROTEIN TECHNOLOGIES INTE...
Tryptophan        Sea lion, Steller, meat with fat (Alaska Native)
Tyrosine         Soy protein isolate, PROTEIN TECHNOLOGIES INTE...
Valine           Soy protein isolate, PROTEIN TECHNOLOGIES INTE...
Name: food, Length: 19, dtype: object

2012년 연방선거관리위원회 데이터베이스

미국연방선거관리위원회는 정치활동 후원금에 대한 데이터를 공개했다. 이 데이터에는 기부자의 이름, 직업, 고용형태, 주소, 기부금액이 포함되어 있다. 그중 재미있는 데이터는 2012년 미국 대통령 선거 데이터다. 2012년 6월에 내려받은 모든 주를 포함하는 전체 제이터는 150메가바이트의 CSV 파일로 파일 이름은 P00000001-ALL.csv이며(이 책의 깃허브 저장소에서 내려받을 수 있다)pandas.read_csv() 함수로 불러올 수 있다.

fec = pd.read_csv('datasets/fec/P00000001-ALL.csv', low_memory=False)

fec.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 1001731 entries, 0 to 1001730
Data columns (total 16 columns):
 #   Column             Non-Null Count    Dtype  
---  ------             --------------    -----  
 0   cmte_id            1001731 non-null  object 
 1   cand_id            1001731 non-null  object 
 2   cand_nm            1001731 non-null  object 
 3   contbr_nm          1001731 non-null  object 
 4   contbr_city        1001712 non-null  object 
 5   contbr_st          1001727 non-null  object 
 6   contbr_zip         1001620 non-null  object 
 7   contbr_employer    988002 non-null   object 
 8   contbr_occupation  993301 non-null   object 
 9   contb_receipt_amt  1001731 non-null  float64
 10  contb_receipt_dt   1001731 non-null  object 
 11  receipt_desc       14166 non-null    object 
 12  memo_cd            92482 non-null    object 
 13  memo_text          97770 non-null    object 
 14  form_tp            1001731 non-null  object 
 15  file_num           1001731 non-null  int64  
dtypes: float64(1), int64(1), object(14)
memory usage: 122.3+ MB

DataFrame에는 다음과 같은 형태로 저장되어 있다.

fec.iloc[123456]
cmte_id             C00431445
cand_id             P80003338
cand_nm         Obama, Barack
contbr_nm         ELLMAN, IRA
contbr_city             TEMPE
                    ...      
receipt_desc              NaN
memo_cd                   NaN
memo_text                 NaN
form_tp                 SA17A
file_num               772372
Name: 123456, Length: 16, dtype: object

기부자와 선거 자금에서 찾을 수 있는 패턴에 대한 통계를 추출하기 위해 이 데이터를 적당한 크기로 쪼개서 나누는 다양한 방법을 떠올릴 수 있을 것이다. 지금까지 배운 내용을 적용한 여러 가지 분석 방법을 살펴보자.

여기에는 정당 가입 여부에 대한 데이터가 없으므로 추가해주는 것이 유용하다. unique() 메서드를 이용해서 모든 정당의 후보 목록을 얻자.

unique_cands = fec.cand_nm.unique()

unique_cands
array(['Bachmann, Michelle', 'Romney, Mitt', 'Obama, Barack',
       "Roemer, Charles E. 'Buddy' III", 'Pawlenty, Timothy',
       'Johnson, Gary Earl', 'Paul, Ron', 'Santorum, Rick',
       'Cain, Herman', 'Gingrich, Newt', 'McCotter, Thaddeus G',
       'Huntsman, Jon', 'Perry, Rick'], dtype=object)
unique_cands[2]
'Obama, Barack'

소속 정당은 dict2를 사용해서 표시할 수 있다.

parties = {'Bachmann, Michelle': 'Republican', 
           'Cain, Herman': 'Republican', 
           'Gingrich, Newt': 'Republican', 
           'Huntsman, Jon': 'Republican',
           'Johnson, Gary Earl': 'Republican',
           'McCotter, Thaddeus G': 'Republican',
           'Obama, Barack': 'Democrat', 
           'Paul, Ron': 'Republican', 
           'Pawlenty, Timothy': 'Republican',
           'Perry, Rick': 'Republican', 
           "Roemer, Charles E. 'Buddy' III": 'Republican',
           'Romney, Mitt': 'Republican', 
           'Santorum, Rick': 'Republican'}

이제 이 사전 정보와 Series 객체의 map() 메서드를 사용해 후보 이름으로부터 정당 배열을 계산해낼 수 있다.

fec.cand_nm[123456:123461]
123456    Obama, Barack
123457    Obama, Barack
123458    Obama, Barack
123459    Obama, Barack
123460    Obama, Barack
Name: cand_nm, dtype: object
fec.cand_nm[123456:123461].map(parties)
123456    Democrat
123457    Democrat
123458    Democrat
123459    Democrat
123460    Democrat
Name: cand_nm, dtype: object
# party 칼럼으로 추가
fec['party'] = fec.cand_nm.map(parties)

fec['party'].value_counts()
Democrat      593746
Republican    407985
Name: party, dtype: int64

분석을 하기 전에 데이터를 다듬어야 한다. 이 데이터에는 기부금액과 환급금액(기부금액이 마이너스인 경우)이 함께 포함되어 있다.

(fec.contb_receipt_amt > 0).value_counts() 
True     991475
False     10256
Name: contb_receipt_amt, dtype: int64

분석을 단순화하기 위해 기부금액이 양수인 데이터만 골라내겠다.

fec = fec[fec.contb_receipt_amt > 0]

버락 오바마와 미트 롬니가 양대 후보이므로 이 두 후보의 기부금액 정보만 따로 추려내겠다.

fec_mrbo = fec[fec.cand_nm.isin(['Obama, Barack', 'Romney, Mitt'])]
fec_mrbo.cand_nm.value_counts()
Obama, Barack    589127
Romney, Mitt     105155
Name: cand_nm, dtype: int64

직업 및 고용주에 따른 기부 통계

직업에 따른 기부 내역 통계는 흔한 조사 방법이다. 예를 들어 변호사는 민주당(Democrat)에 더 많은 돈을 기부하는 경향이 있으며 기업 임원은 공화당(Republican)에 더 많은 돈을 기부하는 경향이 있지만, 그대로 받아들이기보다는 데이터를 통해 직접 확인해보자. 직업별 전체 기부 숫자는 쉽게 구할 수 있다.

fec.contbr_occupation.value_counts()[:10]
RETIRED                                   233990
INFORMATION REQUESTED                      35107
ATTORNEY                                   34286
HOMEMAKER                                  29931
PHYSICIAN                                  23432
INFORMATION REQUESTED PER BEST EFFORTS     21138
ENGINEER                                   14334
TEACHER                                    13990
CONSULTANT                                 13273
PROFESSOR                                  12555
Name: contbr_occupation, dtype: int64

내용을 보면 일반적인 직업 유형이거나 같은 유형이지만 다른 이름으로 많은 결과가 포함되어 있음을 알 수 있다. 아래 코드를 이용해서 하나의 직업을 다른 직업으로 매핑함으로써 이런 몇몇 문제를 제거하자. dict.get()을 사용하는 ‘꼼수’를 써서 매핑 정보가 없는 직업은 그대로 사용한다.

occ_mapping = {
    'INFORMATION REQUESTED PER BEST EFFORTS': 'NOT PROVIDED',
    'INFORMATION REQUESTED': 'NOT PROVIDED', 
    'INFORMATION REQUESTED (BEST EFFORTS)': 'NOT PROVIDED',
    'C.E.O.': 'CEO'
}

# mapping이 없다면 x를 반환한다.
f = lambda x: occ_mapping.get(x, x)
fec.contbr_occupation = fec.contbr_occupation.map(f)

고용주에 대해서도 마찬가지로 처리하자.

emp_mapping = {
    'INFORMATION REQUESTED PER BEST EFFORTS': 'NOT PROVIDED',
    'INFORMATION REQUESTED': 'NOT PROVIDED',
    'SELF': 'SELF-EMPLOYED',
    'SELF EMPLOYED': 'SELF-EMPLOYED',
}

# mapping이 없다면 x를 반환한다.
f = lambda x: emp_mapping.get(x, x)
fec.contbr_occupation = fec.contbr_occupation.map(f)

이제 pivot_table()을 사용해서 정당과 직업별로 데이터를 집계한 다음 최소 2백만불 이상 기부한 직업만 골라내자.

by_occupation = fec.pivot_table('contb_receipt_amt',
                                index='contbr_occupation',
                                columns='party', aggfunc='sum')

over_2mm = by_occupation[by_occupation.sum(1) > 2000000]

over_2mm
party Democrat Republican
contbr_occupation
ATTORNEY 11141982.97 7477194.43
CEO 2074974.79 4211040.52
CONSULTANT 2459912.71 2544725.45
ENGINEER 951525.55 1818373.70
EXECUTIVE 1355161.05 4138850.09
... ... ...
PRESIDENT 1878509.95 4720923.76
PROFESSOR 2165071.08 296702.73
REAL ESTATE 528902.09 1625902.25
RETIRED 25305116.38 23561244.49
SELF-EMPLOYED 741746.40 2245272.64

17 rows × 2 columns

이런 종류의 데이터는 막대그래프(‘barh’는 수평막대그래프를 의미한다)로 시각화하는 편이 보기 좋다(그림 14-12)).

over_2mm.plot(kind='barh')
print('그림 14-12 정당별 최다 기부자들의 직업')
그림 14-12 정당별 최다 기부자들의 직업

오바마 후보와 롬니 후보별로 가장 많은 금액을 기부한 직군을 알아보자. 이 통계를 구하려면 후보 이름으로 그룹을 묶고 이 장의 앞에서 사용했던 변형된 top() 메서드를 사용하면 된다.

def get_top_amouts(group, key, n=5):
    totals = group.groupby(key)['contb_receipt_amt'].sum()
    return totals.nlargest(n)

그리고 직업과 고용주에 따라 집계하면 된다.

grouped = fec_mrbo.groupby('cand_nm')

grouped.apply(get_top_amouts, 'contbr_occupation', n=7)
cand_nm        contbr_occupation    
Obama, Barack  RETIRED                  25305116.38
               ATTORNEY                 11141982.97
               INFORMATION REQUESTED     4866973.96
               HOMEMAKER                 4248875.80
               PHYSICIAN                 3735124.94
                                           ...     
Romney, Mitt   HOMEMAKER                 8147446.22
               ATTORNEY                  5364718.82
               PRESIDENT                 2491244.89
               EXECUTIVE                 2300947.03
               C.E.O.                    1968386.11
Name: contb_receipt_amt, Length: 14, dtype: float64
grouped.apply(get_top_amouts, 'contbr_employer', n=10)
cand_nm        contbr_employer      
Obama, Barack  RETIRED                  22694358.85
               SELF-EMPLOYED            17080985.96
               NOT EMPLOYED              8586308.70
               INFORMATION REQUESTED     5053480.37
               HOMEMAKER                 2605408.54
                                           ...     
Romney, Mitt   CREDIT SUISSE              281150.00
               MORGAN STANLEY             267266.00
               GOLDMAN SACH & CO.         238250.00
               BARCLAYS CAPITAL           162750.00
               H.I.G. CAPITAL             139500.00
Name: contb_receipt_amt, Length: 20, dtype: float64

기부금액

이 데이터를 효과적으로 분석하는 방법은 cut() 함수를 사용해서 기부 규모별로 버킷을 만들어 기부자 수를 분할하는 것이다.

bins = np.array([0, 1, 10, 100, 1000, 10000,
                 100000, 1000000, 10000000])

labels = pd.cut(fec_mrbo.contb_receipt_amt, bins)

labels
411         (10, 100]
412       (100, 1000]
413       (100, 1000]
414         (10, 100]
415         (10, 100]
             ...     
701381      (10, 100]
701382    (100, 1000]
701383        (1, 10]
701384      (10, 100]
701385    (100, 1000]
Name: contb_receipt_amt, Length: 694282, dtype: category
Categories (8, interval[int64, right]): [(0, 1] < (1, 10] < (10, 100] < (100, 1000] < (1000, 10000] < (10000, 100000] < (100000, 1000000] < (1000000, 10000000]]

이제 이 데이터를 이름과 버킷 이름으로 그룹지어 기부금액 규모에 따른 히스토그램을 그릴 수 있다.

grouped = fec_mrbo.groupby(['cand_nm', labels])

grouped.size().unstack(0)
cand_nm Obama, Barack Romney, Mitt
contb_receipt_amt
(0, 1] 493 77
(1, 10] 40070 3681
(10, 100] 372280 31853
(100, 1000] 153991 43357
(1000, 10000] 22284 26186
(10000, 100000] 2 1
(100000, 1000000] 3 0
(1000000, 10000000] 4 0

이 데이터를 보면 오바마는 롬니보다 적은 금액의 기부를 훨씬 많이 받았다. 기부금액을 모두 더한 후 버킷별로 정규화해서 후보별 전체 기부금액 대비 비율을 시각화할 수 있다(그림 14-13).

bucket_sums = grouped.contb_receipt_amt.sum().unstack(0)

normed_sums = bucket_sums.div(bucket_sums.sum(axis=1), axis=0)

normed_sums
cand_nm Obama, Barack Romney, Mitt
contb_receipt_amt
(0, 1] 0.805182 0.194818
(1, 10] 0.918767 0.081233
(10, 100] 0.910769 0.089231
(100, 1000] 0.710176 0.289824
(1000, 10000] 0.447326 0.552674
(10000, 100000] 0.823120 0.176880
(100000, 1000000] 1.000000 0.000000
(1000000, 10000000] 1.000000 0.000000
normed_sums[:-2].plot(kind='barh')
print('그림 14-13 후보별 전체 기부금액 대비 비율')
그림 14-13 후보별 전체 기부금액 대비 비율

기부금액 순에서 가장 큰 2개의 버킷은 개인 후원이 아니므로 그래프에서 제외시켰다.

물론 지금 살펴본 분석은 좀 더 개량할 수 있다. 예를 들어 기부자의 이름과 우편번호를 이용해서 적은 금액을 자주 기부한 사람과 큰 금액을 기부한 사람별로 데이터를 집계할 수도 있을 것이다. 이 책의 저자는 독자들이 이 데이터를 내려받아 직접 살펴보기 강력하게 권장한다.

주별 기부 통계

데이터를 후보자와 주별로 집계하는 것은 흔한 일이다.

grouped = fec_mrbo.groupby(['cand_nm', 'contbr_st'])

totals = grouped.contb_receipt_amt.sum().unstack(0).fillna(0)

totals = totals[totals.sum(1) > 100000]

totals[:10]
cand_nm Obama, Barack Romney, Mitt
contbr_st
AK 281840.15 86204.24
AL 543123.48 527303.51
AR 359247.28 105556.00
AZ 1506476.98 1888436.23
CA 23824984.24 11237636.60
CO 2132429.49 1506714.12
CT 2068291.26 3499475.45
DC 4373538.80 1025137.50
DE 336669.14 82712.00
FL 7318178.58 8338458.81

각 로우를 전체 기부금액으로 나누면 각 후보에 대한 주별 전체 기부금액의 상대적인 비율을 얻을 수 있다.

percent = totals.div(totals.sum(axis=1), axis=0)

percent[:10]
cand_nm Obama, Barack Romney, Mitt
contbr_st
AK 0.765778 0.234222
AL 0.507390 0.492610
AR 0.772902 0.227098
AZ 0.443745 0.556255
CA 0.679498 0.320502
CO 0.585970 0.414030
CT 0.371476 0.628524
DC 0.810113 0.189887
DE 0.802776 0.197224
FL 0.467417 0.532583
percent[:10].plot(kind='barh')
<AxesSubplot:ylabel='contbr_st'>

도시별 기부 통계

데이터를 후보자와 도시별로 집계하였다.

grouped = fec_mrbo.groupby(['cand_nm', 'contbr_city'])

totals = grouped.contb_receipt_amt.sum().unstack(0).fillna(0)

totals[:10]
cand_nm Obama, Barack Romney, Mitt
contbr_city
'WESTPORT 0.0 2500.0
1308 N ASTOR ST C 0.0 500.0
1308 N ASTOR ST CHICAGO 0.0 500.0
29 PALMS 50.0 0.0
33414 STIRRUP LANE 100.0 0.0
4 MARSH RIDGE 200.0 0.0
60154-5908 0.0 75.0
ABBEVILLE 1614.0 4500.0
ABBOTT PARK 0.0 2500.0
ABBOTTSTOWN 145.0 0.0

각 로우를 전체 기부금액으로 나누면 각 후보에 대한 도시별 전체 기부금액의 상대적인 비율을 얻을 수 있다.

percent = totals.div(totals.sum(axis=1), axis=0)

percent[:10]
cand_nm Obama, Barack Romney, Mitt
contbr_city
'WESTPORT 0.000000 1.000000
1308 N ASTOR ST C 0.000000 1.000000
1308 N ASTOR ST CHICAGO 0.000000 1.000000
29 PALMS 1.000000 0.000000
33414 STIRRUP LANE 1.000000 0.000000
4 MARSH RIDGE 1.000000 0.000000
60154-5908 0.000000 1.000000
ABBEVILLE 0.263984 0.736016
ABBOTT PARK 0.000000 1.000000
ABBOTTSTOWN 1.000000 0.000000
percent[:10].plot(kind='barh')
<AxesSubplot:ylabel='contbr_city'>

한국인의 삶 파악하기

대한민국 사람들은 어떻게 살아가고 있을까? 데이터를 분석해 낱낱이 파헤쳐 보자.

‘한국복지패널 데이터’ 분석 준비하기

한국복지패널 데이터는 한국보건사회연구원에서 우리나라 가구의 경제활동을 연구해 복지정책에 반영할 목적으로 발간하는 조사 자료이다. 전국에서 7,000여 가구를 선정해 2006년부터 매년 추적 조사한 자료로, 경제활동, 생활실태, 복지욕구 등 천여 개 변수로 구성되어 있다. 다양한 분야의 연구자와 정책전문가들이 복지패널 데이터를 활용해 논문과 연구보고서를 발표하고 있다.

복지패널 데이터는 엄밀한 절차에 따라 수집되었고 다양한 변수를 담고 있으므로 데이터 분석 기술을 연습하는 데 훌륭한 재료이다. 데이터에는 다양한 삶의 모습이 담겨 있다. 한국복지패널 데이터를 분석해 대한민국 사람들이 어떻게 살아가고 있는지 알아보자.

데이터 분석 준비하기

데이터를 분석하기 위해 준비 작업을 하자.

1. 데이터 준비하기

깃허브에서 Koweps_hpwc14_2019_beta2.sav 파일을 다운로드해 워킹 디렉터리에 삽입한다. 이 파일은 2020년에 발간된 복지패널 데이터로, 6,331가구, 14,418명의 정보를 담고 있다.

2. 패키지 설치 및 로드하기

실습에 사용할 데이터 파일은 통계 분석 소프트웨어인 SPSS 전용 파일이다. pyreadstat 패키지를 설치하면 pandas 패키지의 함수를 이용해 SPSS, SAS, STATA 등 다양한 통계 분석 소프트웨어의 데이터 파일을 불러올 수 있다. 아나콘다 프롬포트에서 다음 코드를 실행해 pyreadstat 패키지를 설치한다.

pip install pyreadstat

노트북을 열고 데이터를 분석하는 데 필요한 패키지를 로드한다.

import pandas as pd
import numpy as np
import seaborn as sns

3. 데이터 불러오기

pd.read_spss()를 이용해 복지패널 데이터를 불러온다. 데이터 원본은 복구할 상황을 대비해 그대로 두고 복사본을 만들어 분석에 활용하겠다.

# 데이터 불러오기
raw_welfare = pd.read_spss('Koweps_hpwc14_2019_beta2.sav')

# 복사본 만들기
welfare = raw_welfare.copy()

4. 데이터 검토하기

데이터의 구조와 특징을 파악해 보겠다.

welfare
h14_id h14_ind h14_sn h14_merkey h_new h14_cobf p14_wsc p14_wsl p14_wgc p14_wgl ... wc14_64 wc14_65 wc14_5aq4 wc14_5aq5 wc14_5aq6 h14_pers_income1 h14_pers_income2 h14_pers_income3 h14_pers_income4 h14_pers_income5
0 2.0 1.0 1.0 20101.0 0.0 NaN 0.291589 0.291589 1307.764781 1307.764781 ... NaN NaN NaN NaN NaN NaN NaN 0.0 NaN
1 3.0 1.0 1.0 30101.0 0.0 NaN 0.419753 0.419753 1882.570960 1882.570960 ... NaN NaN NaN NaN NaN NaN NaN 0.0 NaN
2 4.0 1.0 1.0 40101.0 0.0 NaN 0.265263 0.265980 1189.691668 1192.908537 ... NaN NaN NaN NaN NaN 1284.0 NaN 0.0 NaN
3 6.0 1.0 1.0 60101.0 0.0 NaN 0.494906 0.495941 2219.630833 2224.273816 ... 1.0 . 2.0 4.0 4.0 2304.0 NaN 1800.0 0.0 NaN
4 6.0 1.0 1.0 60101.0 0.0 NaN 1.017935 1.017935 4565.389177 4565.389177 ... 1.0 . 1.0 5.0 2.0 NaN NaN NaN 0.0 NaN
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
14413 9800.0 7.0 1.0 98000701.0 1.0 NaN NaN NaN NaN NaN ... NaN NaN NaN NaN NaN NaN NaN 0.0 NaN
14414 9800.0 7.0 1.0 98000701.0 1.0 NaN NaN NaN NaN NaN ... NaN NaN NaN NaN NaN NaN NaN 0.0 NaN
14415 9800.0 7.0 1.0 98000701.0 1.0 NaN NaN NaN NaN NaN ... NaN NaN NaN NaN NaN 208.0 NaN 0.0 NaN
14416 9800.0 7.0 1.0 98000701.0 1.0 NaN NaN NaN NaN NaN ... 5.0 . 4.0 3.0 3.0 NaN 1200.0 NaN 0.0 NaN
14417 9800.0 7.0 1.0 98000701.0 1.0 NaN NaN NaN NaN NaN ... NaN NaN NaN NaN NaN NaN NaN 0.0 NaN

14418 rows × 830 columns

welfare.shape
(14418, 830)
welfare.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 14418 entries, 0 to 14417
Columns: 830 entries, h14_id to h14_pers_income5
dtypes: float64(826), object(4)
memory usage: 91.3+ MB
welfare.describe()
h14_id h14_ind h14_sn h14_merkey h_new h14_cobf p14_wsc p14_wsl p14_wgc p14_wgl ... wc14_63 wc14_64 wc14_5aq4 wc14_5aq5 wc14_5aq6 h14_pers_income1 h14_pers_income2 h14_pers_income3 h14_pers_income4 h14_pers_income5
count 14418.000000 14418.000000 14418.000000 1.441800e+04 14418.000000 121.000000 11513.000000 11513.000000 11513.000000 11513.000000 ... 2027.000000 2027.000000 2027.000000 2027.000000 2027.000000 2659.000000 3331.000000 989.000000 14418.000000 715.000000
mean 4672.108406 3.121723 1.004855 4.672140e+07 0.201484 2.256198 1.000000 1.000000 4484.952219 4484.952541 ... 3.211643 3.743957 3.513567 4.100641 3.233350 4141.380594 1389.440408 3457.835187 2.038702 1183.292308
std 2792.998128 3.297963 0.143205 2.793014e+07 0.401123 1.675952 0.906021 1.016782 4063.459773 4560.218659 ... 2.174768 3.138629 1.045929 0.937712 1.289456 2583.755449 1211.910836 6619.516319 32.965477 2147.418274
min 2.000000 1.000000 1.000000 2.010100e+04 0.000000 1.000000 0.001998 0.000000 8.960093 0.000000 ... 1.000000 1.000000 1.000000 1.000000 1.000000 0.000000 0.000000 -47000.000000 0.000000 -10600.000000
25% 2356.000000 1.000000 1.000000 2.356030e+07 0.000000 1.000000 0.341814 0.269286 1533.021553 1207.736094 ... 1.000000 2.000000 3.000000 4.000000 2.000000 2448.000000 391.500000 1000.000000 0.000000 206.000000
50% 4535.000000 1.000000 1.000000 4.535010e+07 0.000000 2.000000 0.726304 0.704045 3257.436901 3157.609630 ... 3.000000 3.000000 4.000000 4.000000 3.000000 3540.000000 1116.000000 2498.000000 0.000000 530.000000
75% 6616.000000 7.000000 1.000000 6.616010e+07 0.000000 2.000000 1.366071 1.390045 6126.762919 6234.287538 ... 5.000000 5.000000 4.000000 5.000000 4.000000 5378.500000 2040.000000 4687.000000 0.000000 1295.000000
max 9800.000000 14.000000 9.000000 9.800070e+07 1.000000 7.000000 4.727006 5.790039 21200.393903 25968.049029 ... 9.000000 99.000000 9.000000 9.000000 9.000000 22700.000000 11500.000000 170000.000000 3000.000000 22644.000000

8 rows × 826 columns

지금까지 예제로 사용한 데이터들은 변수의 수가 적고 변수명도 어느정도 이해할 수 있는 단어로 되어 있어서 데이터 구조를 쉽게 파악할 수 있었다. 반면 복지패널 데이터와 같은 대규모 데이터는 변수의 수가 많고 변수명이 코드로 되어 있어 전체 구조를 한눈에 파악하기 어렵다. 규모가 큰 데이터는 데이터 전체를 한 번에 파악하기보다 변수명을 쉬운 단어로 바꾼 다음 분석에 사용할 변수를 하나씩 살펴봐야 한다. 변수를 파악하는 작업은 데이터를 본격적으로 분석하는 과정에서 진행하겠다.

5. 변수명 바꾸기

분석에 사용할 변수 몇 개를 이해하기 쉬운 변수명으로 바꾸겠다. 규모가 큰 조사 자료는 데이터의 특징을 설명해 놓은 코드북codebook을 함께 제공한다. 코드북에는 코드로 된 변수명과 값의 의미가 설명되어 있다. 코드북을 보면 데이터의 특징이 어떠한지 감을 잡을 수 있고, 분석에 어떤 변수를 활용할지, 분석 방향의 아이디어를 얻을 수 있다.

한국복지패널 사이트에서 제공하는 코드북에서 실습에 사용할 변수를 선정해 저자가 Koweps_Codebook_2019.xlsx 파일에 정리해 두었다. 파일을 열어 변수의 특징을 미리 파악하자.

코드북을 참고해 분석에 사용할 변수 7개의 이름을 알아보기 쉬운 단어로 바꾸자.

welfare = welfare.rename(
    columns = {'h14_g3': 'sex',               # 성별
               'h14_g4': 'birth',             # 태어난 연도 
               'h14_g10': 'marriage_type',    # 혼인 상태
               'h14_g11': 'religion',         # 종교
               'p1402_8aq1': 'income',        # 월급
               'h14_eco9': 'code_job',        # 직업 코드
               'h14_reg7': 'code_region'})    # 지역 코드

데이터 분석 절차 살펴보기

데이터를 분석하는데 필요한 준비를 마쳤다. 이제 앞에서 선정한 변수 7개를 이용해 분석하겠다. 여기서는 다양한 분석 주제를 다루는데, 분석마다 두 단계로 진행한다.

1단계 변수 검토 및 전처리

분석에 활용할 변수를 전처리한다. 변수의 특징을 파악하고 이상치와 결측치를 정제한 다음 변수의 값을 다루기 편하게 바꾼다. 전처리는 분석에 활용할 변수 각각 진행한다. 예를 들어 ‘성별에 따른 월급 차이’를 분석한다면 성별, 월급 두 변수를 각각 전처리한다.

2단계 변수 간 관계 분석

전처리를 완료하면 본격적으로 변수 간 관계를 파악하는 분석을 한다. 데이터를 요약한 표와 데이터의 특징을 쉽게 이해할 수 있는 그래프를 만든 다음 분석 결과를 해석한다.

성별에 따른 월급 차이 - 성별에 따라 월급이 다를까?

여성들이 예전에 비해 활발하게 사회 진출을 하고 있지만 직장에서 받는 대우에는 여전히 차별이 있다. 데이터를 분석해 성별에 따라 월급 차이가 있는지 알아보자. 먼저 성별과 월급 두 변수를 검토하고 전처리한 다음 변수 간의 관계를 분석하겠다. 분석 절차를 요약하면 다음과 같다.

분석 절차

1단계: 변수 검토 및 전처리

  • 성별

  • 월급

2단계: 변수 간 관계 분석

  • 성별 월급 평균표 만들기

  • 그래프 만들기

성별 변수 검토 및 전처리하기

1. 변수 검토하기

df.dtypes로 sex(성별) 변수의 타입을 파악하고, df.value_counts()로 범주마다 볓 명이 있는지 알아보자. 출력 결과를 보면 sex는 float64(실수)타입이고 1과 2로 구성된다. 1은 6,505명, 2는 7,913명이 있다.

welfare['sex'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['sex'].value_counts()    # 빈도 구하기
2.0    7913
1.0    6505
Name: sex, dtype: int64

2. 전처리하기

코드북을 보면 성별 변수의 값이 1이면 남자, 2면 여자를 의미한다. 모른다고 답하거나 응답하지 않으면 9로 입력되어 있다. 이 정보를 바탕으로 데이터에 이상치가 있는지 검토하고, 분석할 때 제거하기 편하도록 NaN을 부여해 결측 처리하겠다. 값이 9인 경우도 성별을 알 수 없어 분석에서 제외해야 하므로 결측 처리한다.

# 이상치 확인
welfare['sex'].value_counts()
2.0    7913
1.0    6505
Name: sex, dtype: int64

sex에 1과 2만 있고 9나 다른 값은 없다. 이상치가 없으므로 이상치를 결측 처리하는 절차를 건너뛰어도 된다. 만약 이상치가 있다면 다음과 같이 이상치를 결측 처리하는 작업을 먼저 한 다음, 결측치를 확인해야 한다.

# 이상치 결측 처리
welfare['sex'] = np.where(welfare['sex'] == 9, np.nan, welfare['sex'])

# 결측치 확인
welfare['sex'].isna().sum()
0

sex 변수의 값은 숫자 1과 2로 되어 있다. 값의 의미를 이해하기 쉽도록 문자 ‘male’과 ‘female’로 바꾼 다음 df.value_counts()sns.countplot()를 이요해 바꾼 값이 잘 반영됐는지 출력 결과를 확인하자.

# 성별 항목 이름 부여
welfare['sex'] = np.where(welfare['sex'] == 1, 'male', 'female')

# 빈도 구하기
welfare['sex'].value_counts()
female    7913
male      6505
Name: sex, dtype: int64
# 빈도 막대 그래프 만들기
sns.countplot(data=welfare, x='sex');

성별 변수의 전처리 작업을 완료했다. 이제 같은 절차로 월급 변수를 전처리하자.

월급 변수 검토 및 전처리하기

1. 변수 검토하기

코드북을 보면 월급은 ‘일한 달의 월 평균 임금’을 의미하며 1만 원 단위로 기록되어 있다. income(월급) 면수를 검토하겠다.

성별은 범주 변수이므로 df.value_counts()를 이용해 범주별 빈도를 확인하면 특징을 파악할 수 있었다. 하지만 월급은 연속 변수이므로 df.value_counts()를 이용하면 너무 많은 항목이 출력되어 알아보기 어렵다. 연속 변수는 df.describe()로 요약 통계량을 확인해야 특징을 파악할 수 있다.

welfare['income'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['income'].describe()    # 요약 통계량 구하기
count    4534.000000
mean      268.455007
std       198.021206
min         0.000000
25%       150.000000
50%       220.000000
75%       345.750000
max      1892.000000
Name: income, dtype: float64

출력 결과를 보면 income은 float64 타입이고, 0 ~ 1,892만 원의 값을 지닌다. 150 ~ 345만 원에 가장 많이 분포하고 평균은 268만 원, 중앙값은 평균보다 작은 220만 원으로 전반적으로 낮은 값 쪽으로 치우져 있다.

이번에는 sns.histplot()으로 히스토그램을 만들어 분포를 확인하겠다.

sns.histplot(data=welfare, x='income');    # 히스토그램 만들기

출력된 그래프를 보면 0 ~ 250만 원에 가장 많은 사람이 분포하고, 그 뒤로는 점차 빈도가 감소한다.

2. 전처리하기

코드북을 보면 월급은 만 원 단위로 되어 있고, ‘모름/무응답’은 9999로 코딩되어 있다. 이 정보를 바탕으로 변수를 전처리하겠다.

welfare['income'].describe()    # 이상치 확인
count    4534.000000
mean      268.455007
std       198.021206
min         0.000000
25%       150.000000
50%       220.000000
75%       345.750000
max      1892.000000
Name: income, dtype: float64
welfare['income'].isna().sum()    # 결측치 확인
9884

출력 결과를 보면 최소값이 0, 최대값이 1,892이고, 결측치가 9,884개 있다. 직업이 없어서 월급을 받지 않는 응답자가 있으므로 데이터에 결측치가 있는 것이다. 따라서 월급 변수를 이용해 분석할 때는 먼저 결측치를 제거해야 한다.

코드북에는 ‘모름/무응답’이면 9999로 코딩했다고 되어있다. income의 최대값이 1,892이므로 income이 9999인 행은 없다. 이상치가 없으므로 이상치를 결측 처리하는 절차를 건너뛰어도 된다. 만약 최대값이 9999로 나타나 이상치가 있다면 다음과 같이 이상치를 결측 처리한 다음 결측치가 제대로 만들어졌는지 확인하는 절차를 거쳐야 한다.

# 이상치 결측 처리
welfare['income'] = np.where(welfare['income'] == 9999, np.nan,
                             welfare['income'])

# 결측치 확인
welfare['income'].isna().sum()
9884

성별에 따른 월급 차이 분석하기

1. 성별 월급 평균표 만들기

두 변수의 전처리 작업을 완료했으니 변수 간 관계를 분석할 차례이다. 성별 월급 평균표를 만들어 월급 평균이 성별에 따라 차이가 있는지 비교해 보자.

# 성별 월급 평균표 만들기
sex_income = welfare.dropna(subset=['income']) \
                    .groupby('sex', as_index=False) \
                    .agg(mean_income=('income', 'mean'))

sex_income
sex mean_income
0 female 186.293096
1 male 349.037571

출력 결과를 보면 월급 평균이 남자는 349만 원, 여자는 186만 원으로, 남성이 여성보다 약 163만 원 더 많다.

2. 그래프 만들기

분석 결과를 쉽게 이해할 수 있도록 앞에서 만든 성별 월급 평균표를 이용해 막대 그래프를 만들어보자. 출력된 그래프를 보면 남성의 월급이 여성의 두 배 가까울 정도로 많다.

# 막대 그래프 만들기
sns.barplot(data=sex_income, x='sex', y='mean_income');

나이와 월급의 관계 - 몇 살 때 월급을 가장 많이 받을까?

비정규직이 많아지면서 안정된 직장에 취업하는 것도 어려워졌지만, 젊은 세대를 더욱 힘들게 하는 것은 세대 간 소득 격차가 심해서 사회가 불평등하게 느껴진다는 점이다. 데이터를 분석해서 나이에 따라 월급이 얼마나 다른지 알아보자.

나이 변수를 검토하고 전처리한 다음 나이와 월급의 관계를 분석하겠다. 월급 변수 전처리는 앞에서 완료했으니 생략하겠다.

분석 절차

1단계: 변수 검토 및 전처리

  • 나이

  • 월급

2단계: 변수 간 관계 분석

  • 나이에 따른 월급 평균표 만들기

  • 그래프 만들기

나이 변수 검토 및 전처리하기

1. 변수 검토하기

나이와 월급의 관계를 분석하려면 나이를 나타낸 변수가 있어야 한다. 그런데 한국복지패널 데이터는 나이 변수는 없고 태어난 연도 변수만 있다. 따라서 태어난 연도 변수를 이용해 나이 변수를 만들어야 한다. 먼저 태어난 연도 변수를 검토한 다음 나이 변수를 만들겠다.

welfare['birth'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['birth'].describe()    # 요약 통계량 구하기
count    14418.000000
mean      1969.280205
std         24.402250
min       1907.000000
25%       1948.000000
50%       1968.000000
75%       1990.000000
max       2018.000000
Name: birth, dtype: float64
sns.histplot(data=welfare, x='birth')    # 히스토그램 만들기
<AxesSubplot:xlabel='birth', ylabel='Count'>

2. 전처리

코드북을 보면 태어난 연도는 ‘모름/무응답’일 경우 9999로 코딩되어 있다. 이 정보를 바탕으로 전처리 작업을 하겠다.

welfare['birth'].describe()    # 이상치 확인
count    14418.000000
mean      1969.280205
std         24.402250
min       1907.000000
25%       1948.000000
50%       1968.000000
75%       1990.000000
max       2018.000000
Name: birth, dtype: float64
welfare['birth'].isna().sum()    # 결측치 확인
0

출력 결과를 보면 이상치와 결측치가 없으므로 파생변수를 만드는 단계로 넘어가겠다. 만약 이상치가 발견되면 다음과 같이 전처리한 다음 분석을 진행해야 한다.

# 이상치 결측 처리
welfare['birth'] = np.where(welfare['birth'] == 9999, np.nan,
                            welfare['birth'])

# 결측치 확인
welfare['birth'].isna().sum()
0

3. 파생변수 만들기 - 나이

태어난 연도 변수를 이용해 나이 변수를 만들겠다. 2019년에 조사가 진행됐으니 2019에서 태어난 연도를 뺀 다음 1을 더해 나이를 구하면 된다. 변수를 만들고 df.describe(), sns.histplot()를 이용해 특징을 살펴보겠다.

welfare = welfare.assign(age = 2019 - welfare['birth'] + 1) # 나이 변수 만들기
welfare['age'].describe()
count    14418.000000
mean        50.719795
std         24.402250
min          2.000000
25%         30.000000
50%         52.000000
75%         72.000000
max        113.000000
Name: age, dtype: float64
sns.histplot(data=welfare, x='age');    # 히스토그램 만들기

나이와 월급의 관계 분석하기

월급 변수 전처리는 앞 절에서 완료했다. 나이 변수와 월급 변수의 전처리 작업을 마쳤으니 이제 나이에 따른 월급을 분석할 차례이다.

1. 나이에 따른 월급 평균표 만들기

나이별 월급 평균표를 만들어보자.

# 나이별 월급 평균표 만들기
age_income = welfare.dropna(subset=['income']) \
                    .groupby('age') \
                    .agg(mean_income=('income', 'mean'))

age_income.head()
mean_income
age
19.0 162.000000
20.0 121.333333
21.0 136.400000
22.0 123.666667
23.0 179.676471

2. 그래프 만들기

앞에서 만든 요약표를 이용해 그래프를 만들어보자. x축을 나이, y축을 월급으로 지정하고 나이에 따른 월급의 변화를 나타낸 선 그래프를 만들겠다.

sns.lineplot(data=age_income, x='age', y='mean_income');

출력된 그래프를 보면 20대 초반에 월급을 150만 원가량 받고 이후 지속해서 증가하는 추세를 보인다. 40대에 350만 원가량으로 가장 많이 받고 지속해서 감소하다가 60대 후반부터는 20대보다 낮은 월급을 받는다.

연령대에 따른 월급 차이 - 어떤 연령대의 월급이 가장 많을까?

앞에서는 나이별 월급 평균을 분석했다. 이번에는 나이를 연령대별로 분류한 다음 월급을 비교해 보려고 한다.

분석 절차

1단계: 변수 검토 및 전처리

  • 연령대

  • 월급

2단계: 변수 간 관계 분석

  • 연령대별 월급 평균표 만들기

  • 그래프 만들기

연령대 변수 검토 및 전처리하기

파생변수 만들기 - 연령대

앞에서 만든 나이 변수를 이용해 연령대 변수를 만들겠다. 다음 표의 기준에 따라 연령대 변수를 만든 다음 각 범주에 몇 명이 있는지 살펴보자.

범주 기준
초년층 30세 미만
중년층 30~59세
노년층 60세 이상
# 나이 변수 살펴보기
welfare['age'].head()
0    75.0
1    72.0
2    78.0
3    58.0
4    57.0
Name: age, dtype: float64
# 연령대 변수 만들기
welfare = welfare.assign(ageg=np.where(welfare['age'] < 30, 'young',
                              np.where(welfare['age'] <= 59, 'middle',
                                                             'old')))

# 빈도 구하기
welfare['ageg'].value_counts()
old       5955
middle    4963
young     3500
Name: ageg, dtype: int64
# 빈도 막대 그래프 만들기
sns.countplot(data=welfare, x='ageg');

연령대에 따른 월급 차이 분석하기

월급 변수 전처리는 앞절에서 완료했으니 생략하고 연령대에 따른 월급 차이를 분석하자.

1. 연령대별 월급 평균표 만들기

연령대별로 월급 평균이 다른지 알아보기 위해 연령대별 월급 평균표를 만들어보자. 출력 결과를 보면 월급 평균이 초년층 195만 원, 중년층 329만 원, 노년층 140만 원으로 연령대 별로 차이가 있다.

# 연령대별 월급 평균표 만들기
ageg_income = welfare.dropna(subset=['income']) \
                     .groupby('ageg', as_index=False) \
                     .agg(mean_income=('income', 'mean'))

ageg_income
ageg mean_income
0 middle 329.157157
1 old 140.129003
2 young 195.663424

2. 그래프 만들기

앞에서 만든 요약표를 이용해 그래프를 만들어보자.

# 막대 그래프 만들기
sns.barplot(data=ageg_income, x='ageg', y='mean_income');

막대 정렬 순서는 그래프를 만드는데 사용한 데이터 프레임의 행 순서에 따라 정해진다. 초년, 중년, 노년층 순으로 막대를 정렬하도록 order에 범주 순서를 지정하겠다.

# 막대 정렬하기
sns.barplot(data=ageg_income, x='ageg', y='mean_income',
            order=['young', 'middle', 'old']);

출력된 표와 그래프를 보면 중년층이 330만 원 정도로 가장 많은 월급을 받는다. 노년층의 월급은 140만 원으로, 초년층이 받는 195만 원보다 적다.

연령대 및 성별 월급 차이 - 성별 월급 차이는 연령대별로 다를까?

위에서 성별에 따라 월급 차이가 있는지 분석했었다. 그런데 성별 월급 차이는 연령대에 따라 다른 양상을 보일 수 있다. 이번에는 성별 월급 차이가 연령대에 따라 어떻게 다른지 분석해 보자. 연령대, 성별, 월급 변수 모두 앞에서 전처리 작업을 완료했으므로 전처리 작업은 생략하고 연령대에 따른 성별 월급 차이를 분석해보겠다.

분석 절차

1단계: 변수 검토 및 전처리

  • 연령대

  • 성별

  • 월급

2단계: 변수 간 관계 분석

  • 연령대 및 성별 월급 평균표 만들기

  • 그래프 만들기

연령대 및 성별 월급 차이 분석하기

1. 연령대 및 성별 월급 평균표 만들기

성별 월급 차이가 연령대별로 다른지 알아보기 위해 ‘연령대 및 성별에 따른 월급 평균표’를 만들 것이다.

# 연령대 및 성별 월급 평균표 만들기
sex_income = \
    welfare.dropna(subset=['income']) \
           .groupby(['ageg', 'sex'], as_index=False) \
           .agg(mean_income=('income', 'mean'))

sex_income
ageg sex mean_income
0 middle female 230.481735
1 middle male 409.541228
2 old female 90.228896
3 old male 204.570231
4 young female 189.822222
5 young male 204.909548

2. 그래프 만들기

앞에서 만든 요약표를 이용해 그래프를 만들겠다. 막대가 연령대별로 나열되도록 x축에 ageg를 지정하고, 막대 색깔이 성별에 따라 다르도록 hue에 sex를 지정한다. 축 순서는 order를 이용해 연령대 순으로 설정한다.

# 막대 그래프 만들기
sns.barplot(data=sex_income, x='ageg', y='mean_income', hue='sex',
            order=['young', 'middle', 'old']);

출력된 표와 그래프를 보면 성별 월급 차이의 양상이 연령대별로 다르다. 초년에는 월급 차이가 크지 않다가 중년에 크게 벌어져 남성이 179만 원가량 더 많다. 노년에는 차이가 줄어들지만 여전히 남성이 114만 원가량 더 많다.

앞 절에서 연령대별 월급을 분석할 때 노년층이 초년층보다 월급을 더 적게 받는 것으로 나타났다. 그런데 연령대와 성별로 나눈 이번 분석 결과를 보면 노년층이 초년층보다 월급을 적게 받는 현상은 여성에게서만 타나난다. 남성은 노년층과 초년층의 월급이 비슷하다. 또한, 중년층이 초년층보다 월급을 더 많이 받는 현상은 주로 남성에서 나타나고, 여성은 차이가 크지 않다.

나이 및 성별 월급 차이 분석하기

이번에는 연령대별로 구분하지 않고 ‘나이 및 성별 월급 평균표’를 만들어 선 그래프로 표현하겠다. 그래프에서 성별에 따라 선 색깔이 다르도록 hue에 sex를 입력한다.

# 나이 및 성별 월급 평균표 만들기
sex_age = welfare.dropna(subset=['income']) \
                 .groupby(['age', 'sex'], as_index=False)\
                 .agg(mean_income=('income', 'mean'))

sex_age
age sex mean_income
0 19.0 male 162.000000
1 20.0 female 87.666667
2 20.0 male 155.000000
3 21.0 female 124.000000
4 21.0 male 186.000000
... ... ... ...
140 89.0 male 27.000000
141 90.0 female 27.000000
142 91.0 female 27.000000
143 91.0 male 13.000000
144 92.0 female 27.000000

145 rows × 3 columns

# 선 그래프 만들기
sns.lineplot(data=sex_age, x='age', y='mean_income', hue='sex');

출력된 그래프를 보면 남성의 월급은 50세 전후까지 증가하다가 급격하게 감소하는 반면, 여성은 30세 초반까지 약간 증가하다가 이후로는 완만하게 감소한다. 성별 월급 격차는 30대 중반부터 벌어지다가 50대에 가장 크게 벌어지고, 이후로 점점 줄어들어 80대가 되면 비슷한 수준이 된다.

직업별 월급 차이 - 어떤 직업이 월급을 가장 많이 받을까?

어떤 직업이 월급을 가장 많이 받을까? 직업별 월급을 분석해 보자. 먼저 직업 변수를 검토하고 전처리한다. 월급 변수 전처리 작업은 앞에서 완료했으니 생략하고 직업별 월급 차이를 분석하겠다.

분석 절차

1단계: 변수 검토 및 전처리

  • 직업

  • 월급

2단계: 변수 간 관계 분석

  • 직업별 월급 평균표 만들기

  • 그래프 만들기

직업 변수 검토 및 전처리하기

1. 변수 검토하기

먼저 직업을 나타낸 code_job 변수를 살펴보자.

welfare['code_job'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['code_job'].value_counts()    # 빈도 구하기
611.0    962
941.0    391
521.0    354
312.0    275
873.0    236
        ... 
112.0      2
784.0      2
423.0      1
861.0      1
872.0      1
Name: code_job, Length: 150, dtype: int64

code_job 변수의 값은 직업의 코드를 의미한다. 복지패널 데이터에서 직업은 이름이 아니라 직업분류코드로 입력되어 있다. 지금 상태로는 코드가 어떤 직업을 의미하는지 알 수 없으므로 직업 이름을 나타낸 변수를 만들어야 한다. 한국복지패널 사이트에서 제공하는 코드북의 직업분류코드 목록(한국표준직업분류 제7차 개정)을 보면 코드가 어떤 직업을 의미하는지 알 수 있다.

2. 전처리하기

코드북의 직업분류코드 목록을 이용해 직업 이름을 나타낸 변수를 만들겠다. 먼저 Koweps_Codebook_2019.xlsx 파일의 ‘직종코드’ 시트에 있는 직업분류코드 목록을 불러온 다음 살펴보겠다. 출력 결과를 보면 직업분류코드 목록은 코드와 직업명 두 변수로 구성되고, 직업이 156개로 분류되어 있다.

list_job = pd.read_excel('Koweps_Codebook_2019.xlsx', sheet_name='직종코드')
list_job.head()
code_job job
0 111 의회 의원∙고위 공무원 및 공공단체 임원
1 112 기업 고위 임원
2 121 행정 및 경영 지원 관리자
3 122 마케팅 및 광고∙홍보 관리자
4 131 연구∙교육 및 법률 관련 관리자
list_job.shape  # 행, 열 개수 출력
(156, 2)

df.merge()를 이용해 list_job을 welfare에 결합한다. welfare와 list_job에 공통으로 들어 있는 code_job 변수를 기준으로 결합하면 된다.

# welfare에 list_job 결합하기
welfare = welfare.merge(list_job, how='left', on='code_job')

welfare의 code_job, job 변수 일부를 출력해 잘 결합됐는지 확인해보자. 출력할 때 직업이 결측치인 행은 제외하겠다.

# code_job 결측치 제거하고 code_job, job 출력
welfare.dropna(subset=['code_job'])[['code_job', 'job']].head()
code_job job
2 762.0 전기공
3 855.0 금속기계 부품 조립원
7 941.0 청소원 및 환경미화원
8 999.0 기타 서비스 관련 단순 종사자
14 312.0 경영 관련 사무원

출력 결과를 보면 welfare에 직업 이름으로 된 job 변수가 결합된 것을 확인할 수 있다. 이제 이 변수를 이용해 직업별 월급 차이를 분석하겠다.

직업별 월급 차이 분석하기

월급 변수 전처리는 앞절에서 완료했으니 생략하고 직업별 월급 차이를 분석해보자.

1. 직업별 월급 평균표 만들기

직업별 월급 평균표를 만들겠다. 직업이 없거나 월급을 받지 않는 사람은 분석 대상이 아니므로 제외한다.

# 직업별 월급 평균표 만들기
job_income = welfare.dropna(subset=['job', 'income']) \
                    .groupby('job', as_index=False) \
                    .agg(mean_income=('income', 'mean'))

job_income
job mean_income
0 가사 및 육아 도우미 92.455882
1 간호사 265.219178
2 감정∙기술영업및중개관련종사자 391.000000
3 건물 관리원 및 검표원 168.375000
4 건설 및 광업 단순 종사자 261.975000
... ... ...
142 화학∙고무 및 플라스틱 제품 생산기 조작원 452.714286
143 화학공학 기술자 및 시험원 489.500000
144 환경∙청소 및 경비 관련 관리자 201.000000
145 환경공학∙가스·에너지 기술자 및 시험원 511.000000
146 회계 및 경리 사무원 274.840796

147 rows × 2 columns

2. 그래프 만들기

(1) 월급이 많은 직업

어떤 직업이 월급을 많이 받는지 알아보겠다. 앞에서 만든 요약표를 월급 기준으로 내림차순 정렬하고 상위 10개를 추출한다.

# 상위 10위 추출
top10 = job_income.sort_values('mean_income', ascending=False).head(10)

top10
job mean_income
98 의료 진료 전문가 781.000000
60 법률 전문가 776.333333
140 행정 및 경영 지원 관리자 771.833333
63 보험 및 금융 관리자 734.750000
110 재활용 처리 및 소각로 조작원 688.000000
131 컴퓨터 하드웨어 및 통신공학 전문가 679.444444
24 기계∙로봇공학 기술자 및 시험원 669.166667
6 건설∙전기 및 생산 관련 관리자 603.083333
120 제관원 및 판금원 597.000000
100 의회 의원∙고위 공무원 및 공공단체 임원 580.500000

앞에서 만든 요약표를 이용해 그래프를 만들어보자. 우선 한글로 된 직업 이름이 그래프에 잘 출력되도록 폰트를 설정하겠다.

# 맑은 고딕 폰트 설정
import matplotlib.pyplot as plt
plt.rcParams.update({'font.family': 'Malgun Gothic'})

직업 이름이 길기 때문에 x축에 직업 이름을 지정하면 서로 겹쳐 알아볼 수 없다. 직업 이름을 y축, 월급 평균을 x축에 지정해 그래프를 만들겠다.

# 막대 그래프 만들기
sns.barplot(data=top10, y='job', x='mean_income');

출력된 표와 그래프를 보면 ‘의료 진료 전문가’의 월급이 평균 781만 원으로 가장 많고, 그뒤로는 ‘법률 전문가’, ‘행정 및 경영 지원 관리자’, ‘보험 및 금융 관리자’ 순으로 월급이 많다.

(2) 월급이 적은 직업

이번에는 월급이 적은 직업을 알아보겠다. 직업별 월급 평균표를 월급 기준으로 오름차순 정렬하고 상위 10개를 추출한다.

# 하위 10위 추출
bottom10 = job_income.sort_values('mean_income').head(10)

bottom10
job mean_income
33 기타 돌봄∙보건 및 개인 생활 서비스 종사자 73.964286
34 기타 서비스 관련 단순 종사자 77.789474
128 청소원 및 환경미화원 88.461756
0 가사 및 육아 도우미 92.455882
43 돌봄 및 보건 서비스 종사자 117.162338
97 음식 관련 단순 종사자 118.187500
39 농림∙어업 관련 단순 종사자 122.625000
139 학예사∙사서 및 기록물 관리사 140.000000
126 채굴 및 토목 관련 기능 종사자 140.000000
135 판매 관련 단순 종사자 140.909091

요약표를 이용해 그래프를 만들겠다. 앞에서 만든 월급 상위 10위 그래프와 비교할 수 있도록 그래프의 x축을 0 ~ 800으로 제한하겠다.

# 막대 그래프 만들기
sns.barplot(data=bottom10, y='job', x='mean_income') \
   .set(xlim=(0, 800));

출력된 표와 그래프를 보면 ‘기타 돌봄·보건 및 개인 생활 서비스 종사자’의 월급이 평균 73만 원으로 가장 적고, 그뒤로는 ‘기타 서비스 관련 단순 종사자’, ‘청소원 및 환경미화원’, ‘가사 및 육아 도우미’ 순으로 적다.

월급이 가장 많은 직업과 가장 적은 직업을 비교하면 ‘의료 진료 전문가’는 평균 781만 원, ‘기타 돌봄·보건 및 개인 생활 서비스 종사자’는 평균 73만 원을 받는다. 따라서 ‘의료 진료 전문가’는 ‘기타 돌봄·보건 및 개인 생활 서비스 종사자’의 열 배가 넘는 월급을 받는다.

성별 직업 빈도 - 성별별로 어떤 직업이 가장 많을까?

성 평등이 상식인 세상이지만 여전히 성별에 따라 다른 직업을 선택하는 경향이 있다. 성별에 따라 어떤 직업이 많은지 분석해 보겠다. 성별, 직업 변수 전처리 작업은 앞에서 완료했으니 생략하고 바로 성별과 직업의 관계를 분석하겠다.

분석 절차

1단계: 변수 검토 및 전처리

  • 성별

  • 직업

2단계: 변수 간 관계 분석

  • 성별 직업 빈도표 만들기

  • 그래프 만들기

성별 직업 빈도 분석하기

성별 변수와 직업 변수는 앞절에서 전처리 작업을 완료했으니 생략하고 성별과 직업의 관계를 분석하겠다.

1. 성벽 직업 빈도표 만들기

성별로 직업별 빈도를 구해 상위 10개를 추출하겠다.

# 남성 직업 빈도 상위 10개 추출
job_male = welfare.dropna(subset=['job'])\
                  .query('sex == "male"') \
                  .groupby('job', as_index=False) \
                  .agg(n=('job', 'count')) \
                  .sort_values('n', ascending=False) \
                  .head(10)

job_male
job n
107 작물 재배 종사자 486
104 자동차 운전원 230
11 경영 관련 사무원 216
46 매장 판매 종사자 142
89 영업 종사자 113
127 청소원 및 환경미화원 109
4 건설 및 광업 단순 종사자 96
120 제조 관련 단순 종사자 80
3 건물 관리원 및 검표원 79
141 행정 사무원 74
# 여성 직업 빈도 상위 10개 추출
job_female = welfare.dropna(subset=['job']) \
                    .query('sex == "female"') \
                    .groupby('job', as_index=False) \
                    .agg(n=('job', 'count')) \
                    .sort_values('n', ascending=False) \
                    .head(10)

job_female
job n
83 작물 재배 종사자 476
91 청소원 및 환경미화원 282
33 매장 판매 종사자 212
106 회계 및 경리 사무원 163
31 돌봄 및 보건 서비스 종사자 155
87 제조 관련 단순 종사자 148
73 음식 관련 단순 종사자 126
58 식음료 서비스 종사자 117
88 조리사 114
24 기타 서비스 관련 단순 종사자 97

2. 그래프 만들기

앞에서 만든 직업 빈도표를 이용해 그래프를 만들겠다. 두 그래프를 비교하기 위해 x축 범위를 0 ~ 500으로 통일하겠다.

# 남성 직업 빈도 막대 그래프 만들기
sns.barplot(data=job_male, y='job', x='n').set(xlim=(0, 500));

# 여성 직업 빈도 막대 그래프 만들기
sns.barplot(data=job_female, y='job', x='n').set(xlim=(0, 500));

출력된 표와 그래프를 보면 남녀 모두 ‘작물 재배 종사자’가 가장 많지만, 그뒤로는 순위가 다르다. 남성은 ‘자동차 운전원’, ‘경영 관련 사무원’, ‘매장 판매 종사자’ 순으로 많은 반면 여성은 ‘청소원 및 환경미화원’, ‘매장 판매 종사자’, ‘회계 및 경리 사무원’ 순으로 많다.

종교 유무에 따른 이혼율 - 종규가 있으면 이혼을 덜 할까?

종교가 있으면 이혼을 덜 할까? 종교가 있는 사람이 종교가 없는 사람보다 이혼을 덜 하는지 분석해 보겠다. 먼저 종교, 혼인 상태 두 변수를 검토하고 전처리한 다음 종교와 이혼율의 관계를 분석하겠다.

분석 절차

1단계: 변수 검토 및 전처리

  • 종교

  • 혼인 상태

2단계: 변수 간 관계 분석

  • 종교 유무에 따른 이혼율 표 만들기

  • 그래프 만들기

종교 변수 검토 및 전처리하기

종교 유무를 나타낸 religion을 검토하고 전처리하겠다.

1. 변수 검토하기

welfare['religion'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['religion'].value_counts()    # 빈도 구하기
2.0    7815
1.0    6603
Name: religion, dtype: int64

2. 전처리

앞의 출력 결과를 보면 1과 2 외에 다른 값이 없으므로 이상치를 결측 처리하는 작업은 생략하겠다. 값의 의미를 이해할 수 있도록 코드북의 종교 변수 항목을 참고해 종교 유무를 나타낸 문자를 부여하겠다. 다음 코드의 출력 결과를 보면 종교가 있는 사람이 6,603명, 없는 사람이 7,815명이다.

# 종교 유무 이름 부여
welfare['religion'] = np.where(welfare['religion'] == 1, 'yes', 'no')

# 빈도 구하기
welfare['religion'].value_counts()
no     7815
yes    6603
Name: religion, dtype: int64
# 막대 그래프 만들기
sns.countplot(data=welfare, x='religion');

혼인 상태 변수 검토 및 전처리하기

혼인 상태를 나타낸 marriage_type을 검토하고 전처리하자.

1. 변수 검토하기

welfare['marriage_type'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['marriage_type'].value_counts()    # 빈도 구하기
1.0    7190
5.0    2357
0.0    2121
2.0    1954
3.0     689
4.0      78
6.0      29
Name: marriage_type, dtype: int64

2. 파생변수 만들기 - 이혼 여부

코드북의 ‘혼인 상태’ 변수 항목을 보면 배우자가 있으면 1, 이혼했으면 3으로 입력되어 있다. 이 값을 이용해 이혼 여부를 나타낸 변수를 만들겠다.

내용
0 비해당(18세 미만)
1 유배우
2 사별
3 이혼
4 별거
5 미혼(18세 이상, 미혼모 포함)
6 기타(사망 등)
# 이혼 여부 변수 만들기
welfare['marriage'] = np.where(welfare['marriage_type'] == 1, 'marriage',
                      np.where(welfare['marriage_type'] == 3, 'divorce',
                                                              'etc'))
# 이혼 여부별 빈도
n_divorce = welfare.groupby('marriage', as_index=False) \
                   .agg(n=('marriage', 'count'))

n_divorce
marriage n
0 divorce 689
1 etc 6539
2 marriage 7190
# 막대 그래프 만들기
sns.barplot(data=n_divorce, x='marriage', y='n');

출력한 표와 그래프를 보면 결혼 상태인 사람은 7,190명, 이혼한 사람은 689명이다. 둘 중 어디에도 속하지 않아 ‘etc’로 분류된 사람은 6,539명이다. 이들은 분석 대상이 아니므로 이후 작업에서 제외하겠다.

종규 유무에 따른 이혼율 분석하기

앞에서 전처리한 종교 유무 변수와 이혼 여부 변수를 이용해 분석하겠다.

1. 종교 유무에 따른 이혼율표 만들기

종교 유무에 따른 이혼율표를 만들겠다. 먼저 marriage가 ‘etc’인 경우를 제외하고 ‘종교 유무 및 이혼 여부별 비율’을 구하겠다. value_counts()에 normalize=True를 입력하면 비율을 구한다.

rel_div = welfare.query('marriage != "etc"') \
                 .groupby('religion', as_index=False) \
                 ['marriage'] \
                 .value_counts(normalize=True)

rel_div
religion marriage proportion
0 no marriage 0.905045
1 no divorce 0.094955
2 yes marriage 0.920469
3 yes divorce 0.079531

2. 그래프 만들기

이혼에 해당하는 값만 추출한 다음 proportion을 백분율로 바꾸고 소수점 첫째 자리까지 반올림하겠다. round()는 값을 반올림하는 기능을 한다. round()에 출력할 자릿수를 입력하여 사용한다.

# divorce 추출
# 백분율로 바꾸기
# 반올림
rel_div = rel_div.query('marriage == "divorce"') \
                 .assign(proportion=rel_div['proportion'] * 100) \
                 .round(1)

rel_div
religion marriage proportion
1 no divorce 9.5
3 yes divorce 8.0

앞에서 만든 rel_div를 이용해 막대 그래프를 만들자.

# 막대 그래프 만들기
sns.barplot(data=rel_div, x='religion', y='proportion');

출력한 표와 그래프를 보면 이혼율은 종교가 있으면 8.0%, 종교가 없으면 9.5% 이다. 따라서 종교가 있는 사람이 이혼을 덜 한다고 볼 수 있다.

연령대 및 종교 유무에 따른 이혼율 분석하기

앞에서는 전체를 대상으로 종교 유무에 따른 이혼율을 분석했다. 이번에는 연령대별로 나누어 종교 유무에 따른 이혼율 차이가 연령대별로 어떻게 다른지 알아보자.

1. 연령대별 이혼율표 만들기

우선 이혼율이 연령대별로 다른지 알아보겠다. 연령대 및 이혼 여부별 비율을 구해보자.

age_div = welfare.query('marriage != "etc"') \
                 .groupby('ageg', as_index=False) \
                 ['marriage'] \
                 .value_counts(normalize=True)

age_div
ageg marriage proportion
0 middle marriage 0.910302
1 middle divorce 0.089698
2 old marriage 0.914220
3 old divorce 0.085780
4 young marriage 0.950000
5 young divorce 0.050000

출력 결과를 보면 연령대별로 이혼율이 다르다. 중년층이 8.9%로 가장 높고 그뒤로는 노년층 8.5%, 초년층 5% 순으로 높다. 그런데 빈도를 구해보면 초년층은 결혼, 이혼 모두 다른 연령대에 비해 사례가 적다. 초년층은 사례가 부족해 다른 연령대와 비교하기에 적합하지 않으므로 이후 분석 작업에서 제외하겠다.

# 연령대 및 이혼 여부별 빈도
welfare.query('marriage != "etc"') \
       .groupby('ageg', as_index=False) \
       ['marriage'] \
       .value_counts()
ageg marriage count
0 middle marriage 3552
1 middle divorce 350
2 old marriage 3581
3 old divorce 336
4 young marriage 57
5 young divorce 3

2. 연령대별 이혼율 그래프 만들기

앞에서 만든 age_div에서 초년층을 제외하고 이혼에 해당하는 값만 추출하겠다. 그런다음 proportion을 백분율로 바꾸고 소수점 첫째 자리까지 반올림하겠다.

# 초년 제외, 이혼 추출
# 백분율로 바꾸기
# 반올림
age_div = age_div.query('ageg != "young" & marriage == "divorce"') \
                 .assign(proportion=age_div['proportion'] * 100) \
                 .round(1)

age_div
ageg marriage proportion
1 middle divorce 9.0
3 old divorce 8.6

앞에서 만든 age_div를 이용해 막대 그래프를 만들어보자.

# 막대 그래프 만들기
sns.barplot(data=age_div, x='ageg', y='proportion');

3. 연령대 및 종교 유무에 따른 이혼율표 만들기

종교 유무에 따른 이혼율 차이가 연령대별로 다른지 알아보겠다. 먼저 ‘연령대, 종교 유무, 이혼 여부별 비율’을 구하겠다. 분석에서 초년층은 제외한다.

# etc 제외, 초년 제외
# ageg, religion별 분리
# marriage 추출
# 비율 구하기
age_rel_div = welfare.query('marriage != "etc" & ageg != "young"') \
                     .groupby(['ageg', 'religion'], as_index=False) \
                     ['marriage'] \
                     .value_counts(normalize=True)

age_rel_div
ageg religion marriage proportion
0 middle no marriage 0.904953
1 middle no divorce 0.095047
2 middle yes marriage 0.917520
3 middle yes divorce 0.082480
4 old no marriage 0.904382
5 old no divorce 0.095618
6 old yes marriage 0.922222
7 old yes divorce 0.077778

4. 연령대 및 종교 유무에 따른 이혼율 그래프 만들기

이혼에 해당하는 값만 추출한 다음 proportion을 백분율로 바꾸고 소수점 첫째 자리까지 반올림하겠다.

# divorce 추출
# 백분율로 바꾸기
# 반올림
age_rel_div = \
    age_rel_div.query('marriage == "divorce"') \
               .assign(proportion=age_rel_div['proportion'] * 100) \
               .round(1)

age_rel_div
ageg religion marriage proportion
1 middle no divorce 9.5
3 middle yes divorce 8.2
5 old no divorce 9.6
7 old yes divorce 7.8

앞에서 만든 age_rel_div를 이용해 막대 그래프를 만든다. 종교 유무에 따라 막대 색깔을 다르게 표현하기 위해 hue에 religion을 지정하겠다.

# 막대 그래프 만들기
sns.barplot(data=age_rel_div, x='ageg', y='proportion',
            hue='religion');

출력된 표와 그래프를 보면 중년과 노년 모두 종교가 없는 사람의 이혼율이 더 높다. 중년은 1.3%, 노년은 1.8% 정도로 종교가 없는 사람의 이혼율이 더 높다.

지역별 연령대 비율 - 어느 지역에 노년층이 많을까?

고령 사회가 되면서 노인을 위한 시설을 마련하는 일이 점점 더 중요해지고 있다. 이를 준비하려면 먼저 어느 지역에 노년층이 많이 살고 있는지 알아야 한다. 지역별 연령대 비율을 분석해 어느 지역에 노년층이 많이 사는지 알아보겠다.

먼저 지역 변수를 검토하고 전처리하겠다. 연령대 변수 전처리는 앞에서 완료했으므로 건너뛰고 지역과 연령대의 관계를 분석해보자.

분석 절차

1단계: 변수 검토 및 전처리

  • 지역

  • 연령대

2단계: 변수 간 관계 분석

  • 지역별 연령대 비율표 만들기

  • 그래프 만들기

지역 변수 검토 및 전처리하기

지역을 나타낸 변수 code_region을 검토하고 전처리하자.

1. 변수 검토하기

welfare['code_region'].dtypes    # 변수 타입 출력
dtype('float64')
welfare['code_region'].value_counts()    # 빈도 구하기
2.0    3246
7.0    2466
3.0    2448
1.0    2002
4.0    1728
5.0    1391
6.0    1137
Name: code_region, dtype: int64

2. 전처리하기

code_region 변수의 값은 7개 권역을 의미하는 지역 코드이다. 먼저 코드북을 참고해 지역 코드 목록을 만들겠다. 그런 다음 지역 코드 목록과 welfare에 동시에 들어 있는 code_region 변수를 이용해 welfare에 지역명을 나타낸 변수를 추가하겠다.

내용
1 서울
2 수도권(인천/경기)
3 부산/경남/울산
4 대구/경북
5 대전/충남
6 강원/충북
7 광주/전남/전북/제주도
# 지역 코드 목록 만들기
list_region = pd.DataFrame({'code_region': [1, 2, 3, 4, 5, 6, 7],
                            'region'     : ['서울',
                                            '수도권(인천/경기)',
                                            '부산/경남/울산',
                                            '대구/경북',
                                            '대전/충남',
                                            '강원/충북',
                                            '광주/전남/전북/제주도']})

list_region
code_region region
0 1 서울
1 2 수도권(인천/경기)
2 3 부산/경남/울산
3 4 대구/경북
4 5 대전/충남
5 6 강원/충북
6 7 광주/전남/전북/제주도
# 지역명 변수 추가
welfare = welfare.merge(list_region, how='left', on='code_region')
welfare[['code_region', 'region']].head()
code_region region
0 1.0 서울
1 1.0 서울
2 1.0 서울
3 1.0 서울
4 1.0 서울

지역별 연령대 비율 분석하기

연령대 변수 전처리는 앞에서 완료했으므로 생략하고 지역과 연령대의 관계를 알아보자.

1. 지역별 연령대 비율표 만들기

먼저 ‘지역 및 연령대별 비율표’를 만들겠다.

region_ageg = welfare.groupby('region', as_index=False) \
                     ['ageg'] \
                     .value_counts(normalize=True)

region_ageg
region ageg proportion
0 강원/충북 old 0.459103
1 강원/충북 middle 0.308707
2 강원/충북 young 0.232190
3 광주/전남/전북/제주도 old 0.449311
4 광주/전남/전북/제주도 middle 0.317924
5 광주/전남/전북/제주도 young 0.232766
6 대구/경북 old 0.504051
7 대구/경북 middle 0.296296
8 대구/경북 young 0.199653
9 대전/충남 old 0.413372
10 대전/충남 middle 0.336449
11 대전/충남 young 0.250180
12 부산/경남/울산 old 0.437500
13 부산/경남/울산 middle 0.333742
14 부산/경남/울산 young 0.228758
15 서울 middle 0.385115
16 서울 old 0.376124
17 서울 young 0.238761
18 수도권(인천/경기) middle 0.388170
19 수도권(인천/경기) old 0.325015
20 수도권(인천/경기) young 0.286815

2. 그래프 만들기

region_ageg의 proportion을 백분율로 바꾸고 소수점 첫째 자리까지 반올림하겠다.

# 백분율로 바꾸기
# 반올림
region_ageg = \
    region_ageg.assign(proportion=region_ageg['proportion'] * 100) \
               .round(1)

region_ageg
region ageg proportion
0 강원/충북 old 45.9
1 강원/충북 middle 30.9
2 강원/충북 young 23.2
3 광주/전남/전북/제주도 old 44.9
4 광주/전남/전북/제주도 middle 31.8
5 광주/전남/전북/제주도 young 23.3
6 대구/경북 old 50.4
7 대구/경북 middle 29.6
8 대구/경북 young 20.0
9 대전/충남 old 41.3
10 대전/충남 middle 33.6
11 대전/충남 young 25.0
12 부산/경남/울산 old 43.8
13 부산/경남/울산 middle 33.4
14 부산/경남/울산 young 22.9
15 서울 middle 38.5
16 서울 old 37.6
17 서울 young 23.9
18 수도권(인천/경기) middle 38.8
19 수도권(인천/경기) old 32.5
20 수도권(인천/경기) young 28.7

앞에서 만든 region_ageg를 이용해 지역별 연령대 비율을 나타낸 그래프를 만들어보자. 지역명이 겹치지 않도록 y축에 region, x축에 proportion을 지정해 가로 막대 그래프를 만든다. 연령대별로 막대 색깔을 다르게 표현하도록 hue에 ageg를 지정하자.

# 막대 그래프 만들기
sns.barplot(data=region_ageg, y='region', x='proportion', hue='ageg');

3. 누적 비율 막대 그래프 만들기

앞에서 만든 그래프는 각 지역의 연령대별 비율이 서로 다른 막대로 표현되어 있어서 지역끼리 비교하기 어렵다. 지역끼리 비교하기 쉽도록 연령대별 막대를 누적한 ‘누적 비율 막대 그래프’를 만들겠다.

(1) 피벗하기

행과 열을 회전해 표의 구성을 바꾸는 작업을 피벗pivot이라고 한다. 누적 막대 그래프를 만드는 데 적합하도록 데이터프레임의 행과 열을 회전해 구성을 바꾸겠다.

먼저 그래프를 만드는 데 사용할 지역, 연령대, 비율 변수만 추출한 다음 df.pivot()을 이용해 데이터프레임을 피벗한다.

  • 지역을 기준으로 회전하도록 index=’region’을 입력한다.

  • 연령대별로 열을 구성하도록 columns=’ageg’을 입력한다.

  • 각 항목 값을 비율로 채우도록 values=’proportion’을 입력한다.

다음 코드의 출력 결과를 보면 행은 지역, 열은 연령대로 구성하여 지역 및 연령대별 비율을 나타내고 있다.

# 피벗
pivot_df = \
    region_ageg[['region', 'ageg', 'proportion']].pivot(index='region',
                                                        columns='ageg',
                                                        values='proportion')

pivot_df
ageg middle old young
region
강원/충북 30.9 45.9 23.2
광주/전남/전북/제주도 31.8 44.9 23.3
대구/경북 29.6 50.4 20.0
대전/충남 33.6 41.3 25.0
부산/경남/울산 33.4 43.8 22.9
서울 38.5 37.6 23.9
수도권(인천/경기) 38.8 32.5 28.7

(2) 그래프 만들기

누적 비율 막대 그래프를 만들겠다. 가로 막대 그래프를 만들도록 df.plot.barh()를 이용하고, 막대를 누적하도록 stacked=True를 입력한다.

# 가로 막대 그래프 만들기
pivot_df.plot.barh(stacked=True);

(3) 막대 정렬하기

막대 순서를 정렬하기 위해 pivot_df의 행과 열 순서를 바꾸겠다.

  • 노년층 비율이 높은 순으로 막대 정렬하기: 앞에서 만든 그래프의 막대는 밑에서부터 지역명 가나다순으로 정렬되어 있다. 막대가 노년층 비율이 높은 순으로 정렬되도록 pivot_df의 행을 old 기준으로 정렬하자.

  • 연령대 순으로 막대 색깔 나열하기: 앞에서 만든 그래프는 막대 색깔이 middle(중년), old(노년), young(초년) 순으로 나열되어 있다. 막대 색깔이 초년, 중년, 노년 순으로 나열되도록 pivot_df의 변수 순서를 바꿔보자.

# 노년층 비율 기준 정렬, 변수 순서 바꾸기
reorder_df = pivot_df.sort_values('old')[['young', 'middle', 'old']]

reorder_df
ageg young middle old
region
수도권(인천/경기) 28.7 38.8 32.5
서울 23.9 38.5 37.6
대전/충남 25.0 33.6 41.3
부산/경남/울산 22.9 33.4 43.8
광주/전남/전북/제주도 23.3 31.8 44.9
강원/충북 23.2 30.9 45.9
대구/경북 20.0 29.6 50.4

이제 그래프를 만드는 코드를 실행하면 막대가 노년층 비율이 높은 순으로 정렬되고, 막대 색깔도 연령대 순으로 나열되어 지역끼리 쉽게 비교할 수 있다.

# 누적 가로 막대 그래프 만들기
reorder_df.plot.barh(stacked=True);

출력된 그래프를 보면 ‘대구/경북’의 노년층 비율이 가장 높고, 그 뒤로는 ‘강원/충북’, ‘광주/전남/전북/제주도’, ‘부산/경남/울산’ 순으로 높다.

  1. http://www.ssa.gov/oact/abaynames/limits.html 

  2. 게리 존슨은 나중에 자유당 후보로 지명되긴 했지만 여기서는 공화당 후보로 간주한다. 

댓글남기기