1.请将数据hw1_a和hw1_b分别读入R查看数据并指出各个变量的形式，最小值，最大值，中值，均值，标准差

library(readxl)
hw1_a <- read_excel("Personal/Study/2.R/homework/Practice1/hw1_a.xlsx")
hw1_b <- read_excel("Personal/Study/2.R/homework/Practice1/hw1_b.xlsx")

查看类型：

str(hw1_a)
str(hw1_b)
summary(hw1_a)
summary(hw1_b)

查看最小值，最大值，中值，均值，标准差：

install.packages("plyr")
library(plyr)
each(max,min,median,mean,sd)(form_a$Age)
each(max,min,median,mean,sd)(form_a$Years_at_Employer)
each(max,min,median,mean,sd)(form_a$Years_at_Address)
each(max,min,median,mean,sd)(form_a$Income)

2.将两个数据集进行合并，对于每种数据合并的方式，请说明key并且报告合并后的数据样本总行数

innerjoin<-merge(form_a,form_b)
nrow(innerjoin)
leftjoin<-merge(form_a, form_b, all.x = TRUE)
nrow(leftjoin)
rightjoin<-merge(form_a,form_b,all.y = TRUE)
nrow(rightjoin)
fulljoin<-merge(form_a,form_b,all.x = TRUE,all.y = TRUE)
nrow(fulljoin)

3.请筛选出hw1_a 中收入大于4000的样本，并将此样本和hw1_b 中Is_Default=1的
样本合并，你可以使用inner join的方式。这一问中你可以用pipe的书写形式

hw1_a1=filter(hw1_a,Income>40000)
hw1_b1=filter(hw1_b,Is_Default==1)
inner_join1<-inner_join(hw1_a1,hw1_b1,by="ID")

4.在第2问的基础上, 请给出Income对Years_at_Employer的散点图，你发现了哪些趋
势和现象?

ggplot(data=inner_join)+
  geom_point(mapping = aes(x=Years_at_Employer,y= Income))

5.在第4问的基础上 按照Is_Default 增加一个维度，请展示两变量在不同违约状态
的散点图。请使用明暗程度作为区分方式

ggplot(data=inner_join)+
  geom_point(mapping = aes(x=Years_at_Employer,y= Income,alpha=Is_Default))

6.对于第5问，请使用形状作为另外一种区分方式

ggplot(data=inner_join)+
  geom_point(mapping = aes(x=Years_at_Employer,y= Income,
                           shape=factor(Is_Default)))    

7.请找出各个列的缺失值，并删除相应的行。请报告每一变量的缺失值个数，以及所有缺失值总数

fulljoin[!complete.cases(fulljoin),] #列出存在缺失值的行
newdata<-na.omit(fulljoin)#删除有缺失值的行
newdata
sum(!complete.cases(fulljoin$ID)) #缺失值个数0
sum(!complete.cases(fulljoin$Age)) #缺失值个数11
sum(!complete.cases(fulljoin$Years_at_Employer))#缺失值个数11
sum(!complete.cases(fulljoin$Years_at_Address))#缺失值个数11
sum(!complete.cases(fulljoin$Income))#缺失值个数11
sum(!complete.cases(fulljoin$Credit_Card_Debt))#缺失值个数11
sum(!complete.cases(fulljoin$Automobile_Debt))#缺失值个数11
sum(!complete.cases(fulljoin$Is_Default))#缺失值个数11
sum(is.na(fulljoin))#缺失值总数77个

8.找出Income中的极端值并滤掉对应行的数据

quantile(hw1_a$Income,c(0.025,0.975)) #去掉分位数97.5%以上，2.5%以下的极端值
hw1_a2=filter(hw1_a,Income>14168.81&Income<173030.92)

9.将Income对数化，并画出直方图和density curve.

newdata_outlier$Income <- log(newdata_outlier$Income+0.0001)#收入取对数，真数加一个极小数
ggplot(data = newdata_outlier,mapping = aes(x=newdata_outlier$Income,y=..density..))+
geom_histogram(binwidth=0.1,col="white")+
geom_density(lwd=1,color="orange")

10.以Income作为因变量，Years at Employer作为自变量，进行OLS回归，写出回归的方程，并指出自
变量系数是否在某一显著性水平上显著

m1<-lm(Income~Years_at_Employer,data=hw1_a)
summary(m1)
plot(m1)

回归方程：
Income = 15834 + 3209 * Years_at_Employer
诊断图：
分析：
1.自变量系数和截距P值显著，说明假设：Income和Years_at_Employer线性关系成立，即
Income与Years_at_Employer线性相关
2.R-squared: 0.5584 ，显示模型不够理想，干扰因素较多
3.p-value: < 2.2e-16，说明显著性水平较高

1.编写函数get.root(a,b,c),求解一元二次方程ax²+bx+c=0实根

library(tidyverse)
get.root<-function(a,b,c){
 x1<-NA
 x2<-NA
 if (b^2>4*a*c){
 x1=(-1*b+sqrt(b^2-4*a*c))/(2*a)
 x2=(-1*b-sqrt(b^2-4*a*c))/(2*a)
 print(x1)
 print(x2)
 } else if (b^2==4*a*c)
 {
 x1=(-1*b)/(2*a)
 print(paste("唯一实根",x1))
 } else {
 print("没有实根")
 }}

2.已知某一元二次方程ax²+bx+c=0的三个系数都是随机变量，其中a服从[1，5]上的均匀分布，b服从正态分布N(3,10)，c服从均值1的指数分布。请编写函数get.pob(),计算该方程有实根的概率

a<-runif(100000,min=1,max=5)
b<-rnorm(100000,3,10)
c<-rexp(100000,1/1)
deta<-sum((b^2-4*a*c>=0))
d<-100000
get.pob<-deta/d
get.pob

3.请读入数据，使用软件分别给出 price, marketshare，和brand的缺失值数量。请按照每一个brand, 将数据按照先marjetshare 后price 进行从高到低排序

library(readxl)
data_for_HW2 <- read_excel("Personal/Study/2.R/homework/Practice2/data for HW2.xlsx")
View(data_for_HW2)
sum(is.na(data_for_HW2$price))
sum(is.na(data_for_HW2$marketshare))
sum(is.na(data_for_HW2$brand))

4.请按照brand 的种类，对price和marketshare 求均值

  group_by(brand) %>%
  summarize(pricebar=mean(price))
pricebars
marketsharebars <- data1 %>%
  group_by(brand) %>%
  summarise(marketsharebar=mean(marketshare))
marketsharebars

5.请按照brand 的种类，对price和marketshare 画散点图

  geom_point(mapping = aes(x=price,y=marketshare))+
  facet_wrap(~brand,nrow=2)

6.请按照价格的均值，产生新的变量price_new, 低于均值为“低价格”，高于均值为“高价格”。 同样对市场份额也是，产生变量marketshare_new, 数值为“低市场份额”和“高市场份额”

NAHW$price_new[NAHW$price<=mean(NAHW$price)]<-'低价格'
NAHW$price_new[NAHW$price>mean(NAHW$price)]<-'高价格'
NAHW$marketshare_new[NAHW$marketshare<=mean(NAHW$marketshare)]<-"低市场份额"
NAHW$marketshare_new[NAHW$marketshare>mean(NAHW$marketshare)]<-"高市场份额"

7.请估计模型，marketshare为Y，price为X
8.请画出（7）的拟合直线

m1=lm(marketshare~price,data=data1)
summary(m1)
ggplot(data=data1)+
  geom_point(aes(x=price,y=marketshare))+
  geom_abline(data= m1,col= "blue")

9.请随机产生若干直线，验证（5）的结果是最优的

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]<-(marketshare[i]-b0[n]-b1[n]*price[i])^2}
  sum[n]<-sum(d)
  n<-n+1
}

resi=m1$residuals
resi2=sum(resi^2)

check=sum(as.numeric(sum<resi2))

10.请估计模型，marketshare为Y，price和brand 为X

LM1<-lm(marketshare~price+brand,NAHW)
LM1
summary(LM1)
#P值0.4589>0.05表明不够显著
ggplot(NAHW, aes(price, marketshare)) +
 geom_point(aes(color = brand)) +
 geom_smooth(se = FALSE)

1.A 和 B 约定在某篮球场见面。他俩都不太守时，出现时间服从均匀分布。他俩也都没有耐心， 每个人都会只等对方十分钟就会离开。已知 A 到篮球场的时间为下午 4 点到 5 点之间。
(1) 如果 B 到达篮球场的时间也为下午 4 点到 5 点之间，模拟运行 50000 次，看看他们成功相遇的概率。
(2) 对上一问的 50000 次模拟，用不同颜色在一张图中展示成功相遇与否。
(3) B 应该如何选择 4 点到 5 点之间的哪个时间段，来提升他们成功相遇的概率？ 用模拟展示你的理由
（1)模拟运行 50000 次，AB成功相遇的概率

 meet<-function(n){m=0
 for(i in 1:n) {x<-round(runif(1,0,3600))
                    y<-round(runif(1,0,3600))
              if(abs(x-y)<=600){m=m+1}}
print(m/n)}
meet(50000)

A与B各在3600秒中取一个时间点，如果两人的时间差的绝对值小于等于600秒，则为一次相遇，标记一次，然后循环50000次，其中相遇的次数占总次数的比例则为相遇概率
(2)用不同颜色在一张图中展示成功相遇与否

x<-round(runif(50000,0,3600))
y<-round(runif(50000,0,3600))
z<-abs(x-y)<=600
library(tidyverse) 
data<-tibble(x,y,z)    
ggplot(data=data)+
  geom_point(mapping = aes(x=x,y=y,color=z))

（3）B 应该如何选择 4 点到 5 点之间的哪个时间段，来提升他们成功相遇的概率？

p<-function(y){if(y<600){p<-(600+y)/3600}

  else if(y>3000){p<-(3600-(y-600))/3600}

  else{p<-1200/3600}

  return(p)}

p(500)

p(3200)

p(700)

如第二问绘图所示，相遇蓝色区域可分为三段

x<600

600≤x≤3000

x>3000

分别对各段求相遇概率

以上可知在600≤x≤3000（即A在4点10分-4点50分）他们成功相遇的概率最高
**2.请使用 nycflights13 和 pipe 语法
（1）从 flights 数据表中挑选出以下变量：(year, month, day, hour, origin, dep_delay, distance, carrier)，将生产的新表保存为 flight1

library(dplyr)
library(nycflights13)
flight1<-flights%>%
 select(year,month,day,hour,origin,dep_delay,distance,carrier) 

（2）从 weather 数据表中挑选出以下变量： (year, month, day, hour, origin, humid, wind_speed)，将生产的新表保存为 weather1

weather1<-weather%>%
  select(year,month,day,hour,origin,humid,wind_speed)

（3）将 flight1 表和 weather1 表根据共同变量进行内连接，随机抽取 100000 行数据，将生产的结果保存为 flight_weather。 (提示：sample_n()函数，不用重复抽取)

    flight_weather <- inner_join(flight1, weather1) %>% sample_n(100000)

（4）从 flight_weather 表中对三个出发机场按照平均出发延误时间排降序，并将结果保留在 longest_delay 表中。把结果展示出来。

longest_delay<-flight_weather %>%
 group_by(origin) %>%
  summarise(meanDELAY=mean(dep_delay,na.rm = TRUE))

（5）根据出发地（origin） 在同一个图中画出风速 wind_speed（x 轴）和出发延误时间dep_delay（y 轴） 的平滑曲线图。

ggplot(data = flight_weather)+
geom_smooth(mapping = aes(x = wind_speed, y = dep_delay, linetype = origin))

（6）根据不同出发地（origin）在平行的 3 个图中画出风速 wind_speed（x 轴）和出发延误时间 dep_delay（y 轴）的散点图

ggplot(data = flight_weather) + 
geom_point(mapping = aes(x = wind_speed, y = dep_delay)) + facet_wrap(~origin,nrow = 1)

（7）根据 flight_weather 表，画出每个月航班数的直方分布图，x 轴为月份，y 轴是每个月份航班数所占的比例。

ggplot(data=flight_weather)+geom_bar(mapping = aes(month, y=..prop.., group = 1) )

（8）根据 flight_weather 表，画出每个月航班距离的 boxplot 图，x 轴为月份，y 轴为航行距离, 根据的航行距离的中位数从低到高对 x 轴的月份进行重新排序。

ggplot(data=flight_weather)+geom_boxplot(mapping = aes(x=reorder(month, distance, FUN=median), y=distance))

image.png

(1)

(H <- function(p) -sum(p*log(p)))

(2)

(DKL <- function(p,q) sum( p*(log(p)-log(q)) ))

(3)

IB <- list()
IB[[1]] <- c( 0.2 , 0.2 , 0.2 , 0.2 , 0.2 )
IB[[2]] <- c( 0.8 , 0.1 , 0.05 , 0.025 , 0.025 )
IB[[3]] <- c( 0.05 , 0.15 , 0.7 , 0.05 , 0.05 )
purrr::map_dbl( IB , H )
[1] 1.6094379 0.7430039 0.9836003

(4)

Dm <- matrix( NA , nrow=3 , ncol=3 )
for ( i in 1:3 ) {
  for ( j in 1:3 ) {
    Dm[i,j] <- DKL( IB[[j]] , IB[[i]] )
  }
} 
Dm