Author Archives: ys2887

Python learning form

Capture.JPGCapture.JPGCapture.JPG

Capture.JPG

Capture.JPG

Capture.JPG

parameter can not be changed but arguments can

parameter is used in a def function while argument is used in the existing function.

Capture.JPG

 

Capture.JPGCapture.JPGCapture.JPG

Capture.JPG

Capture.JPGCapture.JPG

list

Capture.JPG

while and for

Capture.JPG

the use of function

Capture.JPG

Fibonacci number

Capture.JPG

sort function

Capture.JPG

 

 

Capture

CaptureCapture.JPG

 

deal with indentation in terminal:

notice that when the loop ends there should be a blank row by an enter

Screen Shot 2016-08-19 at 12.02.46 PM.png

Screen Shot 2016-08-19 at 12.02.54 PM.png

Matlab tips

R tips

在通常的情况下，我习惯写一些比较有意思的感悟，这样喜欢看的受众也广泛一些。但这次恐怕要让部分朋友失望了，我想总结下统计软件R的学习精要。R软件做为使用最为广泛的开源统计软件，有着其独特的设计和吸引人的魅力。但这次与大家谈的，更多为R软件的入门精要，期望对统计感兴趣或者工作中用到统计的朋友能够快速上手，体会到用R来完成数据分析工作的乐趣。

整篇文章分为两个部分。在基础部分会向大家介绍R的一些基本概念，包括：R是什么、相关资源和安装、数据结构、实用函数、输入输出、作图、脚本编程。而在实战部分会向大家介绍如何用R来进行一些常用的数据分析工作：相关系数、回归分析、分类判别、聚类分析、时间序列、主成分分析、T校验。

基础部分 

  1. R是什么？英文的解释为：The R Project for Statistical Computing，它既是一种统计语言，也是一种统计软件。
  2. 为什么用R？总结了如下四点原因：

a)  多领域的统计资源：目前在R网站上约有2400个程序包，涵盖了基础统计学、社会学、经济学、生态学、空间分析、系统发育分析、生物信息学等诸多方面。

b)  跨平台：R可在多种操作系统下运行，如Windows、MacOS、多种Linux和UNIX等。

c)  命令行 和 脚本编程:R即时解释，输入命令，即可获得相应的结果。也支持相对复杂功能的脚本编程，处理输入输出文件等。

d)   开源，庞大的社区支持:可以免费、持续的获取最新的算法支持。

3. R的资源

a)         官方网站：http://www.r-project.org/

b)         中文网站：http://www.r-user.org/

c)         统计之都： http://cos.name/

4. 安装流程：共分为四步

a)         下载源代码版本（Source Code for all Platforms）

b)         配置安装路径：./configure –prefix=/home/work/tools/R/

c)         编译：make

d)         安装：make install

5. 数据结构：R与常见的编程语言一样，具有数值、字符、复数、逻辑四种基本类型。同时，为了方便的应对统计和分析处理，提供了矩阵、数据框等复合类型来表示二维表格。常用的复合类型如下，可以用函数mode()查看一个变量的类型。

a)         向量：一个变量

b)         因子：一个分类变量

c)         数组：你懂的

d)         矩阵：相当于2维数组

e)         数据框：相当于含有多种数据类型元素的矩阵，而矩阵元素必须是同一类型

f)          时间序列：附加时间、频率等额外的数据

g)         列表：相当于结构体，里面什么都可以有。R很多函数的返回值都是这个类型，可以用  列表变量名$内部变量名，来取出该变量内容。

6.  R的一些实用函数

a)         ？ lm：查看该函数的说明文档，用于已知函数名时的精准查询

b)         apropos(“lm”)： 模糊查询所有相关函数

c)         attributes(x)：查看列表变量x所含有的变量名，相应变量结果可以用‘x$变量名’得到。通常用来查看函数返回结果中包含哪些可用内容。

d)         summary(x)：查看列表变量x的汇总内容，通常用于查看函数的返回结果

e)         setwd()、getwd()：设置工作路径

f)          library(MASS)：加载一些非默认的函数库。R将所有的函数库分为三个层次，最为通用+常用的自动加载到内存，可以直接在程序中调用。而次要不太常用的放在硬盘上，需要使用的时候需要用library函数现加载到内存中。最专用、只能满足小众需求的，比如DNA生物库，是不在默认安装包里的，需要使用时可以用install.packages()来下载安装。

7. 输入输出：最常用的为读入输出变量，和读入输出表格。

a)         直接写变量名：打印到屏幕

b)         从文件输入

i.              X=scan(“file.txt”,nrow=3,ncol=5,byrow=T)  #将文件内容读入变量

ii.              X=read.table(“file.txt”, header=T, sep=“\t”)  #读入表格，数据有表头，以”\t”作为列分隔符

iii.              X=read.delim(“file.txt”, header=F, sep=”\t”)  #读入表格，但每行的列数可以不同

c)         写入文件

i.              write(X，file=“output”,append=F)

ii.              write.table(X，file=“output”，append=T，row.names=F，col.names=F，sep=”\t”)

8.   作图：R将作图分为高级作图函数和低级作图函数，其中高级函数提供的是标准、成型的图形服务，如散点图、柱状图等。而低级函数提供基本的绘画功能，高级作图函数是调用低级函数来实现的。如果喜欢，也可以用低级函数编写自己想要的图形。

9.   最后再进入实战之前，介绍下写程序所需的最后一部分知识，即程序的控制和执行。R软件的条件语句和循环控制语句与C语言一致，使用很方便。

a)         R的条件判断：If（condition） statement1  else  statment2

b)         R的循环语句：for（i  in  1:n） { statement }   while（condition） { statement }

c)         保存脚本文件为 program.R，执行脚本： R  -f  “program.R”

 

实战部分

上一部分为R的基本使用知识，但只有这些我们并没有必要去学习R，任何一种编程语言均能实现这些基本的运算，在大部分方面甚至可以实现的更好。那么R吸引我们的魅力就要在这一部分向大家展现，即实用的数学、统计、算法函数。

1.       回归分析

回归是寻找自变量和因变量之间的逻辑关系的方法。R提供的常用回归分析有：线性回归lm()、非线性回归nls()、Logistic回归glm()。线性回归和非线性回归，需要给出一个设计好的公式，而Logistic回归则要求回归变量为分类变量。下面以一个实例来阐述R软件的处理流程：

a)         读入数据文件，为每行有6列数据的表格

b)         对每一行数据进行线性回归，取得结果中的斜率项

c)         将这行数字、数字之和、线性拟合的数据输出到结果文件中

#共六个数据项

x<-c(1:6)

#读入数据文件

data = read.delim(“input”,header=FALSE,sep=”\t”)

for(i in 1:length(data[,1]))

{

#取出一列数据

y<-as.numeric(data[i,])

cl<-data.frame(x,y)

#最简单的线性公式

 result<-lm(y~x, data=cl)

#复杂自设计的非线性公式

#result1<-nls(y~a+(0.49-a)*exp(-b*(x-8)), data=cl, start=list(a=0.1, b=0.01))

#Logistic回归，回归结果变量为定性变量

#result2=glm(Opinion~Age+Income, family=”binomial”, data=cl)

#取出结果中的线性斜率的系数项

rate=as.matrix(coef(result))

#打印结果，附加求和项和系数项

result=c(data[i,], sum(data[i,]),rate[2,1])

write.table(result,file=”output”,  append=TRUE,row.names=FALSE,col.names=FALSE,sep=”\t”)

}

 

2.       相关系数

线性相关系数为确定两个变量是否数据相关的产用方法，R实现了三种算法：pearson、kendall、spearman。输入文件为两列，第一列的title是first，第二列的title是second。

x=read.table(“input”, head=T)

#三种线性相关系数均支持, 默认为Pearson系数

result=cor(x$first, x$second)

#result=cor(x$first, x$second,method=”kendall”)

#result=cor(x$first, x$second,method=”spearman”)

 

Result

 

3.       聚类分析

聚类是自动根据数据特征进行分类的方法，其特点时事先并不知道分类依据，完全根据自变量数据表现分类。常用方法除了大家都知道的KMeans聚类之外，还有分层聚类：将每步都将聚类最近的点进行合并，直到合并成只有一个节点，根据需要取某一步为最终结果。

x=read.table(“cities1.txt”, head=T, sep=”,”)

#得到数量为4的聚类结果

result=kmeans(x[,-1], 4, iter.max = 10)

attributes(result)

result$cluster

result$centers

#分层聚类的用法，计算距离的方案有: 离差平方和”ward”, 最短距离法”single”, 最长距离法”complete”, 类平均法”average”, “mcquitty”, 中间距离法”median” or 重心法”centroid”.

hc=hclust(dist(x[,-1]),”ave”)

#在linux下执行，会将画图结果保存到当前目录下，为PDF文件

plot(hc, label=x[,1])

 

4.       判别分析

分类的方案有很多，这里仅仅取一个用线性模型判别的方法。需要注意，分类与聚类不同，需要事先给出一个训练集，同时需要预留出测试集来对分类效果进行测试。

#该数据分析函数不在R默认包中，需要引入

library(MASS)

#数据集使用内置数据集 iris

iris

#设定随机种子

set.seed(5)

#从150个数据记录中取一半作为训练集，一半作为测试集

trainset = c(sample(1:50, 25), sample(51:100,25), sample(101:150, 25))

#使用线性判别模型进行判别, 将Species的列作为判别结果，其余列作为输入

model = lda(Species~., data=iris, subset=trainset)

#使用训练的模型对测试集进行判别

result=predict(model, iris[-trainset,])

#打印判别的结果

result$class

 

5.       时间序列分析

对时间序列分析有两个常用的目标：第一个是拆解，将数据中的趋势变量，周期变量拆解出来，以便对数据的趋势有更清晰的认识。第二个是预测，预估曲线未来的发展趋势。

R中的曲线拆解将曲线分为了三个部分：趋势分量、周期分量、剩余部分。而我的个人观点，将任何一条时间序列曲线拆成如下四个部分趋势+季节性+例外（异常）+随机波动更加合理一些。另：ts函数是将普通的数字串组织成时间序列的格式。

x=scan(“input”)

#转换为时间序列的格式

#start参数或者为数字，或者为年月的二元组

charge=ts(x,frequency=7,start=1)

#时间序列的分解函数

a=stl(charge, “period”)

result=a$time.series

#将分解结果的三部分，输出到文件中  第一列为季节分量，第二列为趋势分量，第三列为剩下的部分

write.table(result,file=”output”, append=FALSE,row.names=FALSE,col.names=FALSE,sep=”\t”)

#ts.plot(a$time.series)

 

曲线趋势预估，最常用的方法有移动平均法、指数平滑法、ARIMA方法，其中移动平均法非常简单，这里仅提供R软件关于指数平滑法和ARIMA方法的使用程序。

x=scan(“input”)

charge=ts(x,frequency=7,start=1)

b = HoltWinter(charge, beta=0)

#输出拟合的曲线

b$fitted

#对未来三十天作预测

predict(b, n.ahead=30)

 

ARIMA预测模型，需要自行根据自相关图和偏相关图来确定P与Q，函数会自动训练系数值。

x=scan(“input”)

charge=ts(x,frequency=7,start=1)

#查看曲线的自相关图和偏自相关图

#acf(charge)   pacf(charge)

#arima模型拟合

#第二项参数为p,d,q 分别为自回归模型项数,差分等级,移动平均项数,seasonal的三个参数为季节性分量的配置

w=arima(charge, c(0,1,1),seasonal=list(order=(1,2,1), period=7)

#模型的拟合系数

w$coef

#使用模型预测

result=predict(w, n.ahead=30)

write.table(result,file=”output”, append=FALSE,row.names=FALSE,col.names=FALSE,sep=”\t”)

 

 

6.       主成分分析

确定一个多维数据的主成分，将多维数据进行降维又能最大化的保留原始信息。

#使用R软件自带的iris数据集

iris

#进行主成分分析

#算法使用相关矩阵还是协方差矩阵来计算

prin_result=princomp(iris[c(2,3,4)], cor=T)

prin_result$sdev

prin_result$loadings

#因子分析，只选择一个因子

factor_result=factanal(iris[c(2,3,4)], fact=1)

#载荷结果

factor_result$loadings

 

主成分分析结果，有两部分信息是最为重要的。第一是每个新的主成分所占原数据信息的比重，如下面comp.1的比重是1.49是最大的；comp.2的比重是0.86，也还可以；而最后comp.3则没有多大分量，可以忽略。在这个数据结果下，我们可以取前两个合成变量做为数据的主成分。Loading载荷表示每个合成变量如何由原始变量来生成，如comp.1=0.418*Sepal.width+petal.length*-0.648+petal.width*-0.636。

 Comp.1    Comp.2    Comp.3

1.4904641 0.8624409 0.1863127

Loadings:

Comp.1 Comp.2 Comp.3

Sepal.Width   0.418  0.907

Petal.Length -0.648  0.256 -0.717

Petal.Width  -0.636  0.335  0.695

 

 

7.       T校验

最后给大家介绍下统计中最简单最基本的T校验（即通过统计实验来判断一个结论是否是真的），让大家轻松下吧~

#读入数据，厂商声明每袋糖的平均重量是500克，实际观测的50袋重量倾向于小于500克

x=scan(“sugar.txt”)

#假设整体数据分布为正态分布，且中心为500，得到这样观测的概率p

#由于实际观测的数据倾向于小于500克，所以用“less”的单边校验

result=t.test(x, m=500,alt=”less”)

#对两个样本总体均值差异的校验，可用于验证产品试验效果

#result=t.test(data1，data2)

result$p.value

#如果p小于1%，即说明这种情况是小概率事件，有理由认为厂商的声明是错误的

if (result$p.value<0.01)

{

“the assumption is wrong!”

}

 

最近我面试了好多有意向来百度的同学，我深切的感觉到：学校中学习知识的思路太刻板了，同学们都绞尽脑汁的去回忆各种算法或概念的公式。这恰恰是忽略了事物的本质，学习一个知识或算法，就要站在这个知识或算法的发明者的角度去思考，究竟是为了解决怎样的问题？而这个方案中运用怎样的原理解决问题的？或许这就是工作后的学习更加务实了吧，所以期望这篇日志能对日常用统计来处理数据的同学有一定的帮助。第二部分的程序都是直接可运行来查看结果的，有需要相关实验数据的可以单独联系我。