Hmisc 包中的 describe()函数 ______包的_______函数可返回变量和观测的数量、info 值、缺失值和唯一值的数目、平均值、分位数，以及五个最大的值和五个最小的值
reshape 包/reshape2 包 ________包是一套重构和整合数据集的万能工具
vcd 包 Arthritis 数据集；
sm 包 sm.density.compare(x, factor) _______包中的_______函数可向图形叠加两组或更多的核密度图
pastecs 包中有一个名为 stat.desc()的函数
Hmisc 包中的 describe()函数 ______包的_______函数可返回变量和观测的数量、info 值、缺失值和唯一值的数目、平均值、分位数，以及五个最大的值和五个最小的值
psych 包的 describe()的函数 ______包的_______函数可以计算非缺失值的数量、平均数、标准差、中位数、截尾均值、绝对中位差、最小值、最大值、值域、偏度、峰度和平均值的标准误。
doBy 包 的summaryBy doBy 包中的_______函数可以分组计算概述统计量
psych 包 describe.by
vcd 包的 assocstats()函数 _______包中的_______函数可以用来计算二维列联表的 phi 系数、列联系数和 Cramer’s V 系数
ggm 包的 pcor()函数 _______包的________函数可以计算偏相关系数
psych 包的 corr.test() _______包的________函数可以为 Pearson、Spearman 或 Kendall 相关计算相关矩阵，同时报告显著性水平
原文链接：http://tecdat.cn/?p=13033
介绍
布丰投针是几何概率领域中最古老的问题之一。它最早是在1777年提出的。它将针头掷到有平行线的纸上，并确定针和其中一条平行线相交的可能性。令人惊讶的结果是概率与pi的值直接相关。
R程序将根据上段所述的情况估算pi的值并使用gganimate进行动态可视化。
第1部分
对于A部分，我们创建一个数据帧，该数据帧将在3个不同的区间上生成随机值，这些区间将代表x，y的范围以及每个落针点的角度。这是一个易于实现的随机数情况，需要使用runif函数。此功能要求输入数量，后跟一个间隔。生成数字后，我们会将值保存到数据框中。
rneedle <- function(n) {
x = runif(n, 0, 5)
y = runif(n,0, 1)
angle = runif(n,-pi, pi) #从-180到180的角度
values<-data.frame(cbind(x, y, angle))
return(values)
}
values<-rneedle(50)
#检查是否生成50×3矩阵
values
#我们的数据帧已经成功生成。
x y angle
1 4.45796267 0.312440618 1.3718465
2 3.43869230 0.462824677 2.9738367
3 2.55561523 0.596722445 -2.9638285
4 3.68098572 0.670877506 -0.6860502
5 0.03690118 0.202724803 -0.3315141
6 4.64979938 0.180091416 -0.3293093
7 4.92459238 0.172328845 -0.5221133
8 3.50660347 0.752147374 2.9100221
9 2.03787919 0.167897415 -0.3213833
10 0.38647133 0.539615776 -0.1188982
11 3.28149935 0.102886770 -1.6318256
12 3.68811892 0.765077533 1.2459037
13 1.52004894 0.682455494 -0.4219802
14 3.76151379 0.508555610 0.1082087
...
第2部分
我们绘制第一部分中的针。重要的是不要在这个问题上出现超过2条水平线。它使我们可以进行检查以了解此处描绘的几何特性的一般概念。话虽如此，让我们注意我们决定在每个方向上将图形扩展1个单位。原因是想象一个针尾从y = 1开始，其角度为pi / 2。我们需要假设该方向的范围最大为2。
plotneedle(values)
第3部分
在下面，将基于阅读布冯针和基本几何原理的知识，查看pi的估算值。
buffon(values)
第4部分
运行代码后，我们收到以下答案。
>buffon(X)
[1] 3.846154
set.seed(10312013)
X <- rneedle(50)
plotneedle(X)
buffon(X)
>buffon(X)
[1] 3.846154
第5部分
如前几节所述，当我们投掷更多的针头时，我们期望以最小的不确定性获得更准确的答案。从Approxpi函数运行代码后，我们收到了平均值= 3.172314和方差0.04751391的值。对于这样一个简单的实验，它对pi进行了很高的估计。
Approxpi(500)
mean(Approxpi(500))
var(Approxpi(500))
>mean(Approxpi(500))
[1] 3.172314
>var(Approxpi(500))
[1] 0.04751391
接下来对模拟次数从500~600的预测进行动态可视化,红色表示针投放到了直线上：
参考资料
Schroeder，L.（1974年）。布冯针问题：许多数学概念的激动人心的应用。
最受欢迎的见解
1.R语言动态图可视化：如何、创建具有精美动画的图
2.R语言生存分析可视化分析
3.Python数据可视化-seaborn Iris鸢尾花数据
4.r语言对布丰投针（蒲丰投针）实验进行模拟和动态
5.R语言生存分析数据分析可视化案例
6.r语言数据可视化分析案例：探索brfss数据数据分析
7.R语言动态可视化：制作历史全球平均温度的累积动态折线图动画gif视频图
8.R语言高维数据的主成分pca、 t-SNE算法降维与可视化分析案例报告
9.python主题LDA建模和t-SNE可视化
R语言-v1-基础知识
Iretara 12-17 21:18
以例题的形式简述R语言基础知识
# 读取文件
setwd(" 文件链接的时候，用 / ")
install.packages(" readxl ")
library(readxl)
library (tidyverse)
hw1_a<- read_excel ("hw1_a.xlsx", col_types=c("numeric", "numeric", "numeric", "numeric", "numeric") )
hw1_b<- read_excel ("hw1_b.xlsx")
#读取csv
library(readr)
hw1_a<- read_csv ("/")
View(hw1_a)
# 描述型函数
hw1_a + hw1_b 表
#描述最小值，最大值，中值，均值，标准差
Str (hw1_a) #查看数据并指出各个 变量的形式
summary (hw1_a) #指出各个变量的形式， 最小值，最大值，中值，均值
library(psych)
describe (hw1_a) #比summary更简便的方法， 可以直接读取标准差等；但是，使用describe不可读取 NA值, 可以尝试使用 Hmisc包中 describe
描述型函数-R
# 连接
hw1_a %>% inner_join (hw1_b, by ="ID")
hw1_a %>% left_join (hw1_b, by ="ID")
hw1_a %>% right_join (hw1_b, by ="ID")
hw1_a %>% full_join (hw1_b, by ="ID")
inner_join<- inner_join (hw1_a,hw1_b, by =“ID”) #报告合并后的 总行数 ，178行
full_join<- full_join (hw1_a,hw1_b, by ="ID")
( nrow (full_join)) #报告合并后的 总行数 ，200行
> length (full_join$ID)
#找出各个列的 缺失值
i<-NA
a<-NA
for(i in 1:length(full_join[1,])){ a[i]<- sum(is.na( full_join[,i] ) ) }
paste("缺失值是",a)
#缺失值总数
sum(is.na(full_join))
#删除缺失值 na.omit()
full_join1=filter(full_join,!is.na(full_join[2]))
full_join1=filter(full_join1,!is.na(full_join1[3]))
full_join1=filter(full_join1,!is.na(full_join1[4]))
full_join1=filter(full_join1,!is.na(full_join1[5]))
full_join1=filter(full_join1,!is.na(full_join1[6]))
full_join1=filter(full_join1,!is.na(full_join1[7]))
full_join1=filter(full_join1,!is.na(full_join1[8]))
sum(is.na(full_join1))
找出Income中的 极端值 并滤掉对应行的数据
quantile (hw1_a$Income,c(0.025,0.975))
hw1_a2= filter (hw1_a,Income>14168.81 &Income<173030.92)
#使用dplyr进行数据转换
arrange()
>arrange (hw1_a,Income) #默认升序
>arrange(hw1_a, desc (Income)) #desc降序，NA排序一般最后
select()
>select (hw1_a, - (Years_at_Address:Income)) #不要变量
>rename (hw1_a, In_come=Income) #改名
>select（hw1_a,Income, exerything ()） #把Income放在前面
拓例题1：
library(nycflights13)
view(flights)
#counts
（1）
not_cancelled <- flights %>%
filter(! is.na(dep_delay), !is.na(arr_delay))
（2）
not_cancelled %>%
group_by (year,month,day) %>%
summarize (mean=mean(dep_delay))
（3）
delays <- not_cancelled %>%
group_by (tailnum) %>%
summarize (delay=mean(arr_delay))
ggplot (data=delays,mapping=aes(x= delay))+
geom_freqpoly (binwidth=10) #freqpoly
（4）
delays <- not_cancelled %>%
group_by(tailnum) %>%
summarize(delay=mean(arr_delay,na.rm=TRUE), n=n() ) #tailnum的次数
ggplot(data=delays,mapping=aes(x= n, y=delay))+
geom_point(alpha=1/10)
拓例题2：
#请按照价格的均值，产生新的变量price_new, 低于均值为“低价格”，高于均值为“高价格”。 同样对市场份额也是，产生变量marketshare_new, 数值为“低市场份额”和“高市场份额”
price=data1$price
pricebar=mean(price)
price_new= ifelse (price>pricebar,“高价格”,”低价格”)
marketshare=data1$marketshare
marketsharebar=mean(marketshare)
marketshare_new=ifelse(marketshare>marketsharebar ,“高市场份额”,”低市场份额”)
data1= mutate (data1,price_new,marketshare_new)
#可视化
#将Income 对数化
lninc<- log (hw1_a$Income)
#画出直方图和 density curve密度曲线
hist (lninc,prob=T)
lines ( density (lninc),col="blue")
# 添加额外变量 的办法，在 aes()中添加 样式 (color、size、alpha、shape)
ggplot(data=inner_join)+
geom_point(mapping = aes(x=Years_at_Employer,y= Income, alpha= Is_Default))
# 按照Is_Default 增加一个维度，使用明暗程度作为区分方式
ggplot(data=inner_join)+
geom_point(mapping = aes(x=Years_at_Employer,y= Income,
alpha=factor( Is_Default ) ))
#使用形状作为另外一种区分方式
ggplot(data=inner_join)+
geom_point(mapping = aes(x=Years_at_Employer,y= Income,
shape=factor( Is_Default)))
可视化-R
拓展：
#将 flight1 表和 weather1 表根据共同变量进行内连接，随机抽取 100000 行数据， 将生产的结果保存为 flight_weather。 (提示：sample_n()函数，不用重复抽取)
flight_weather <- inner_join(flight1, weather1) %>% sample_n(100000)
# 从 flight_weather表中对三个出发机场按照平均出发延误时间排降序，并将结果保留在 longest_delay表中。把结果展示出来
longest_delay<- flight_weather %>%
group_by(origin) %>%
summarize(delay=mean(dep_delay, na.rm=TRUE )) %>%
arrange(desc(delay))
#根据不同出发地（origin）在平行的 3 个图中画出风速 wind_speed（x 轴）和出发 延误时间 dep_delay（y 轴）的散点图。
ggplot(data= flight_weather) +
geom_point(mapping=aes(x=wind_speed,y=dep_delay))+
facet_grid(.~origin, nrow = 3 ) # 按照class分类，分成3行
#根据 flight_weather 表，画出每个月航班数的直方分布图，x 轴为月份，y 轴是每个 月份航班数所占的比例。
ggplot(data=flight_weather)+
geom_bar(mapping=aes(x=month, y=..prop .., group=1))
#根据 flight_weather 表，画出每个月航班距离的 boxplot 图，x 轴为月份，y 轴为 航行距离, 根据的航行距离的中位数从低到高对 x 轴的月份进行重新排序
ggplot(data=flight_weather)+
geom_boxplot(mapping=aes(x= reorder (month,distance,FUN=median),y=distance))
线性回归
# 以Income作为因变量，Years at Employer作为自变量，进行 OLS回归
m1<- lm (Income ~ Years_at_Employer,data=hw1_a)
#通过***判断显著性
summary (m1)
#画出拟合直线
ggplot(data= hw1_a)+
geom_point(aes(x=Income,y=Years_at_Employer))+
geom_abline(data= m1,col= "blue")
#证明拟合直线是最优的
b0=runif(20000,-5,5)
b1=runif(20000,-5,5)
d<-NA
sum<-NA
n<-1
while(n<=20000){
for(i in 1:24){
d[i]<-(hw1_a $ Income[i]-b0[n]-b1[n]*hw2$ Years_at_Employer[i])^2}
sum[n]<-sum(d)
n<-n+1
}
resi=m1$residuals
resi2=sum(resi^2)
check=sum(as.numeric(sum<resi2))
check
