R로 된 feature selection 하는 코드입니다.
어디선가 퍼왔는데 출처가 기억이 나질 않네요.
부연설명을 드리면 모델링을 할 때 feature(또는 독립변수)중 어떤 것이 중요한 것인지 판별하는 코드입니다.
http://github.com/euriion/code_snippets/blob/master/R/feature_selection.R
카테고리 보관물: 컴퓨터언어
Python multi core 구동 코드
Python을 이용해서 ETL의 일부인 파싱이나 전처리 작업을 수행하는 경우가 많습니다.
빅데이터인 경우에도 데이터를 Hadoop이나 Hive 또는 Oracle과 같은 RDBMS에 로딩하기 전에 할 수 있는 것들은 최대한 전처리를 한 후에 사용하는 경우가 많이 있습니다.
물론 데이터량이 아주 많으면 Map/Reduce를 작성하는 것이 더 낫습니다만 그리 크지 않은 데이터는 한 대의 서버에서 자원을 풀가동해서 처리해 버리는 것이 작업속도를 줄일 수 있습니다.
Hadoop이 일반화되기 이전에는 이런 형태의 코드를 더 구체화해서 여러 대의 서버에서 동시에 구동되도록 (마치 맵리듀스처럼) 프로세스를 돌리고 결과를 취합하는 것을 만드는 것이 빈번했었습니다.
https://gist.github.com/euriion/5719443
코드를 수정하면 더 복잡한 것도 할 수 있습니다만 매우 복잡하다면 다른 구조를 생각해 보는 것이 좋습니다.
CSV포맷을 TSV포맷으로 바꾸는 간단한 스크립트
엑셀(Excel)에서 CSV 포맷으로 파일을 저장할 때 텍스트 컬럼을 Escaping처리하는 경우가 있습니다.
주로 쉼표(comma)와 따옴표(double quotation)을 그렇게 변환해 버리는데 Hadoop이나 이 포팻을 Hive에 업로드해서 사용하려면 Escaping을 빼야 합니다.
크기가 크지 않은 CSV는 간단하게 Python으로 변환코드를 작성해서 올려서 사용하는 것이 편한데 그럴때 사용했던 소스코드입니다.
R에서 데이터를 로딩할 때도 이 방법이 편합니다.
이런 간단한 작업도 넓은 의미에서는 데이터 먼징 (Data Munging) 포함됩니다.
https://gist.github.com/euriion/5720809
R ARIMA 예제 코드
R의 ARIMA 모형의 예제입니다.
서버의 메모리의 사용량의 추이를 보고 얼마 후에 고갈되는지를 예측하는 코드입니다.
물론 예측력은 많이 떨어지고 현실성이 없을 수 있습니다.
# -------------------------
# Memory usage forecasting
# -------------------------
library(stats)
arima(lh, order = c(1,0,0))
arima(lh, order = c(3,0,0))
arima(lh, order = c(1,0,1))
arima(lh, order = c(3,0,0), method = "CSS")
arima(USAccDeaths, order = c(0,1,1), seasonal = list(order=c(0,1,1)))
arima(USAccDeaths, order = c(0,1,1), seasonal = list(order=c(0,1,1)),
method = "CSS") # drops first 13 observations.
# for a model with as few years as this, we want full ML
arima(LakeHuron, order = c(2,0,0), xreg = time(LakeHuron)-1920)
## presidents contains NAs
## graphs in example(acf) suggest order 1 or 3
require(graphics)
(fit1 <- arima(presidents, c(1, 0, 0)))
tsdiag(fit1)
(fit3 <- arima(presidents, c(3, 0, 0))) # smaller AIC
tsdiag(fit3)
# ----- prediction part
od <- options(digits=5) # avoid too much spurious accuracy
predict(arima(lh, order = c(3,0,0)), n.ahead = 12)
(fit <- arima(USAccDeaths, order = c(0,1,1),
seasonal = list(order=c(0,1,1))))
predict(fit, n.ahead = 6)
options(od)
# ----- Arima
library(forecast)
fit <- Arima(WWWusage,c(3,1,0))
plot(forecast(fit))
x <- fracdiff.sim( 100, ma = -.4, d = .3)$series
fit <- arfima(x)
plot(forecast(fit,h=30))
# ----- Arima forecast for memory usage (unit %) -----
library(forecast) # need to install the package "forecast"
memory.usage.threshold <- 100 # 100%
memory.usage.forecast.period <- 30 # 미래 30일분까지 예측
memory.usage.observations.startdate <- "2012-09-01"
memory.usage.observations <- c(10,11,30,35,36,39,48,56,75,69,68,72,71,72,83) # 관측치 12일분
memory.usage.period <- seq(as.Date(memory.usage.observations.startdate), length=length(memory.usage.observations), by="1 day") # 날짜세팅
memory.usage.df <- data.frame(row.names=memory.usage.period, memory=memory.usage.observations) # data.frame으로 변환
memory.usage.ts <- ts(data=memory.usage.df) # time series 생성
memory.usage.model <- auto.arima(memory.usage.ts) # arima 모델 생성
memory.usage.forecast <- forecast(memory.usage.model, h=memory.usage.forecast.period) # forecast 결과 생성
memory.usage.forecast.df <- as.data.frame(memory.usage.forecast) # forecast 결과 변환
d = memory.usage.threshold,][1,])) # 100 이 넘는 최초 데이터 추출
if(is.na(d)) {
print(sprintf("앞으로 %s일동안 %s%% 초과하지 않음", memory.usage.forecast.period, d - length(memory.usage.observations)))
} else {
print(sprintf("%s일 후에 %s%% 초과됨", d - length(memory.usage.observations), memory.usage.threshold))
}
# ---- 시각화(Plotting)
plot(memory.usage.forecast) # plotting
abline(h=100, col = "red", lty=3)
abline(v=d, col = "red", lty=3)
library(ggplot2)
library(scales)
plt <- ggplot(data=pd,aes(x=date,y=observed))
p1a<-p1a+geom_line(col='red')
p1a<-p1a+geom_line(aes(y=fitted),col='blue')
p1a<-p1a+geom_line(aes(y=forecast))+geom_ribbon(aes(ymin=lo95,ymax=hi95),alpha=.25)
p1a<-p1a+scale_x_date(name='',breaks='1 year',minor_breaks='1 month',labels=date_format("%b-%y"),expand=c(0,0))
p1a<-p1a+scale_y_continuous(name='Units of Y')
p1a<-p1a+opts(axis.text.x=theme_text(size=10),title='Arima Fit to Simulated Datan (black=forecast, blue=fitted, red=data, shadow=95% conf. interval)')
원본 소스코드는 아래에 있습니다.
https://github.com/euriion/code_snippets/blob/master/R/forecast_exam.R
R ggplot2 – 경제인구동향 그래프 찍기
R에서 ggplot2로 경제활동인구찍기를 해봤습니다.
사실은 다른 것을 플로팅해보려다가 원하는 자료를 다운로드 받는 것이 만만치 않아서
대충 지나가다가 통계청 데이터중에 처음 보는 것을 가지고 찍어본 것입니다.
우선은 데이터를 가져와야 합니다. 통계청에 가시면 여러가지 통계데이터를 제공하고 있습니다.
아래 사이트에 가서 경제활동인구동향데이터를 긁어 옵니다.
http://kosis.kr/feature/feature_0103List.jsp?mode=getList&menuId=03&NUM=180
CSV로 다운로드 받아서 해도 되겠지만 데이터가 크지 않으므로 그냥 소스코드에 집어넣기 위해서 copy&paste를 해버립니다.
사이트에서 바로 복사하면 컬럼간의 구분이 Tab으로 되어 있을텐데요.
편집기에서 제가 Tab문자를 쓰지 않아서 Tab을 모두 세미콜론(;)으로 바꿨습니다. 그리고 header를 month와 population으로 해서 column 이름을 아예 데이터에서 지정해버렸습니다.
코드는 아래와 같습니다.
economic_activity_population <- "month;population
2009.09;24,630
2009.10;24,655
2009.11;24,625
2009.12;24,063
2010.01;24,082
2010.02;24,035
2010.03;24,382
2010.04;24,858
2010.05;25,099
2010.06;25,158
2010.07;25,232
2010.08;24,836
2010.09;24,911
2010.10;25,004
2010.11;24,847
2010.12;24,538
2011.01;24,114
2011.02;24,431
2011.03;24,918
2011.04;25,240
2011.05;25,480
2011.06;25,592
2011.07;25,473
2011.08;25,257
2011.09;25,076
2011.10;25,409
2011.11;25,318
2011.12;24,880
2012.01;24,585
2012.02;24,825
2012.03;25,210
2012.04;25,653
2012.05;25,939
2012.06;25,939
2012.07;25,901
2012.08;25,623"
자 이제 data.frame으로 로딩합니다.
statdata <- read.table(file=textConnection(econonic_activity_population), header = TRUE, sep = ";", quote = ""'", as.is=TRUE,colClasses=c("character", "character"))
그리고 나서
statdata$month 컬럼을 Date형식으로 바꿔줍니다. 그러지 않으면 나중에 곤란해집니다. 궁금하시면 직접해 보시구요.
년도와 날짜로만 되어 있는 문자열을 날짜형으로 바꾸기 위해서 강제로 01을 붙여서 그달의 첫째날로 바꿔버립니다.
그리고 바꿀때 타임존(tz)을 서울(Asia/Seoul)로 해줍니다. 안해주면 가끔 날짜가 UTC로 바뀌는 경우가 있습니다.
statdata$month <- as.Date(paste(statdata$month, ".01", sep=""), "%Y.%m.%d", tz="Asia/Seoul")
그리고 statdata$population을 숫자형으로 바꿔줍니다. 그런데 숫자데이터에 콤머가 있으므로 콤머를 다 제거해주고 숫자형으로 바꿉니다.
코드는 아래와 같습니다.
statdata$population <- as.numeric(gsub(",", "", statdata$population))
자 이제 플로팅을 해보죠. ggplot2를 로딩한 다음에 바로 찍습니다. ggplot2를 설치하지 않으셨으면 먼저 설치하셔야 합니다.
install.packages("ggplot2") # 설치를 안했으면 먼저 설치부터...
library(ggplot2) # 로딩
ggplot(statdata, aes(x=month, y=population)) + geom_line()
플로팅한 그림은 아래와 같습니다.
나왔네요. 그런데 회색배경에 검은선이라 이쁘지 않네요.
선에 색을 넣어 봅니다.
ggplot(statdata, aes(x=month, y=population)) + geom_line(colour="blue")
조금 낫긴하지만 역시 허전합니다.
아래가 비어서 그런것 같으니 선을 그리지 말고 채우기를 이용해서 그려보죠.
사실 이런류의 데이터는 색을 채우지 않고 선으로 보는 것이 데이터를 보기에는 더 좋습니다만 누군가에게 던져 줄때는 예쁘게 출력하는 것도 중요합니다.
본인이 보는 그래프는 선으로 된 것이면 충분하지만 누군가에게 보여주는 그래프는 알록달록해야 합니다.
gplot(statdata, aes(x=month, y=population)) + geom_area()
찍었습니다. 그런데 어라? 그래프의 모양이 바뀐것 같습니다. 자세히 보니 Y축이 영역(range)가 바뀐것 같은데 그것 때문에 밋밋해진 것이군요.
geom_area는 기본적으로 Y축의 0값부터 채워 넣기 때문에 geom_line과는 다르게 반응합니다.
그래서 Y축의 레인지를 강제로 고쳐줘야 합니다. 밑부분의 넙다란 부분을 잘라서 없애버리는 것입니다.
ggplot(statdata, aes(x=month, y=population)) + geom_area() + coord_cartesian(ylim = c(23500, 26500))
그럴듯하게 나오네요. 그런데 여전히 색이 단조롭네요.
테두리의 선색과 채울때 쓴 색을 다르게 줘서 입체감을 조금 살려보죠. 조금 옅은 회색과 조금 진한 회색을 써보겠습니다.
ggplot(statdata, aes(x=month, y=population)) + geom_area(colour="gray10", fill="gray50") + coord_cartesian(ylim = c(23500, 26500))
조금 괜찮아졌습니다만 역시 흑백보다는 칼라를 넣는 것이 좋을 것 같네요.
삷이 우울해질것만 같습니다.
제가 좋아하는 보라색을 넣어봅니다.
ggplot(statdata, aes(x=month, y=population)) + geom_area(colour="gray10", fill="gray50") + coord_cartesian(ylim = c(23500, 26500))
유후~ 괜찮아졌네요.
자. 이제 대충 색은 되었고 나머지를 조금 더 손을 봐보죠.
Y축의 값은 단위를 천명으로 한 숫자입니다. 통계청 데이터의 설명에 그렇게 나와 있습니다.
원래 숫자대로 바꿔보죠. 값에 곱하기 1000을 해서 원래 숫자로 바꿉니다.
귀찮으니 이쯤에서 처음부터 데이터를 다시 로딩해서 바꿔버립니다. 기존 데이터에 그냥해도 무방합니다.
statdata <- read.table(file=textConnection(econonic_activity_population), header = TRUE, sep = ";", quote = ""'", as.is=TRUE,colClasses=c("character", "character"))
statdata$month <- as.Date(paste(statdata$month, ".01", sep=""), "%Y.%m.%d", tz="Asia/Seoul")
statdata$population <- as.numeric(gsub(",", "", statdata$population)) * 1000
네 바꿨습니다. Y축 레인지를 조절하는 값도 1000을 곱해서 바꿔줘야 하는것을 기억하시구요.
이제 Y축의 레이블값의 숫자들에게 콤머를 찍어줍니다. 콤머를 찍으려면 scales 패키지를 로딩해야 합니다.
library(scales)
ggplot(statdata, aes(x=month, y=population)) + geom_area(colour="#5c0ab9", fill="#8a4fcd") + coord_cartesian(ylim = c(23500000, 26500000)) + scale_y_continuous(labels=comma)
이제 X축의 레이블들을 수정해봅니다. 2011-01 처럼 되어 있는 것을 2011년 011월 이렇게 바꿔줄 것입니다.
한글 폰트를 지정해 주시는 것을 잊지 말아야 합니다.. 저는 맥을 사용하므로 애플산돌고딕네오를 적었습니다만 Windows 사용자라면 “맑은 고딕”같은 것을 사용해주세요. 폰트를 정확히 지정하지 않으면 한글이 출력되지 않습니다.
ggplot(statdata, aes(x=month, y=population)) + geom_area(colour="#5c0ab9", fill="#8a4fcd") + coord_cartesian(ylim = c(23500000, 26500000)) + scale_y_continuous(labels=comma) + scale_x_date(labels = date_format("%Y년 %m월")) + theme(axis.text.x = element_text(family="Apple SD Gothic Neo"))
자! 되었습니다.
이제 마지막으로 X축과 Y축의 타이틀을 한글로 바꿔줍니다.
month와 population을 년/월과 경제활동인구라는 말로 바꿔줄 것입니다.
그러면서 코드를 좀 깔끔하게 정리해 보죠. 한줄에 너무 덕지덕지 다 붙여 놓으면 나중에 찾아서 고치기가 어렵습니다.
ggp <- ggplot(statdata, aes(x=month, y=population))
ggp <- ggp + geom_area(colour="#5c0ab9", fill="#8a4fcd")
ggp <- ggp + coord_cartesian(ylim = c(23500000, 26500000))
ggp <- ggp + scale_y_continuous(labels=comma)
ggp <- ggp + scale_x_date(labels = date_format("%Y년 %m월"))
ggp <- ggp + xlab("년도/월") + ylab("경제활동인구")
theme.title <- element_text(family="Apple SD Gothic Neo", face="bold", size=12, angle=00, hjust=0.54, vjust=0.5)
theme.text <- element_text(family="Apple SD Gothic Neo", size=10)
ggp <- ggp + theme(axis.title.x = theme.title, axis.title.y = theme.title, axis.text.x = theme.text)
ggp
rm(theme.title)
rm(theme.text)
gc()
드디어 완성입니다.
아주 이쁘지는 않지만 그럭저럭 보여줄만한 수준이 되는 것 같습니다.
다음 포스트에서는 플롯의 영역별로 색을 지정하는 것을 해보죠.