R數據科學###

20200919 8:30 raining####

shif+enter 換行

planes %>%
count(tailnum)%>%
filter(n>1)

weather%>%
count(year,month,day,hour,origin)%>%

后面的年月日等條件捆綁計數

filter(n>1)

flights2%>%
select(-origin,-dest)%>%
left_join(airlines,by="carrier")

left_join左連接

x%>%
right_join(y,by="key")

right_join右連接

x%>%
inner_join(y,by="key")

inner_join內連接 取交集

x%>%
full_join(y,by="key")

full_join全鏈接 取并集

semi_join(x,y,by="key")

semi_join保留x表中與Y表中的觀測相匹配的所有x

anti_join(x,y,by="key")

anti_join去除x表中與Y表中的觀測相匹配的所有x

top_dest <-flights%>%
count(dest,sort=TRUE)%>%

sort=TRUE降序

head(10)

head()前幾位

top_dest

優(yōu)化的方法optim()

best<-optim(c(0,0),measure_distance,data=sim1)
best$par

sim1_mod<-lm(y~x,data=sim1)

lm()線性模型

coef(sim1_mod)

coef()調出

線性回歸

普通最小二乘

20200920 8:30 經管院b224 raining####

y是均值為0滥崩，標準差為2計算10000次的正態(tài)分布栗柒，x是男女 隨機取10000次帮辟，y=b+ax

y<-rnorm(10000, mean=0, sd=2)
x<-sample(
0:1,10000,replace=TRUE
)
mm<-lm(y~x)

mm<-lm(y~x-1) 減掉1表示去除截距

coef(mm)

textA=scan("D:/mikistudysoftware/RStudio/text A.txt",what="",sep="",na.strings="NA")

導入文本文件馍悟，注意文本路徑""改成"/"

what:聲明讀入為字符類型數據评姨，可能指定讀入的精度/類型

sep：指定各個讀入的數據之間的分隔符帮坚；默認情況下分隔符：空格盏道、tab稍浆；如果不是其它分隔符，例如“：/”通過SEP來指定

可以通過list指定讀入變量的變量名猜嘱，同時生成的對象為列表衅枫，則可以同時讀入字符與數字

Skip 從第幾行開始讀入數據

Nlines 指定最大讀入行數

如果通過鍵盤輸入的時候，不typeof(hw1_a$ID)希望出現下標提示朗伶，則可以使用：quiet=TRUE

encoding =””指定的編碼格式弦撩，有時候讀入的中文可能會出現亂碼的時候，可能通過這個參數來指定：Latin-1 或者 UTF-8

juhaoA=sum(str_detect(textA,"\."))

統計文本A中的.的個數

wordA=length(textA)

統計文本A文字的個數

willA=sum(str_detect(textA,"will"))

統計文本A中“will”這個詞語的個數

20200927

箱線圖適用于連續(xù)變量異常值檢測

箱線圖

ggplot(hw1_fj_dna,aes(,Income))+

geom_boxplot()

去除極端值的箱線圖

boxplot(hw1_fj_dna$Income, outline = F)

boxplot(hw1_fj_dna$Income)

極端值

boxplot.stats(hw1_fj_dnaout

去除極端值

hw1_fj_dna_IncomeFilter<-hw1_fj_dnaIncome%in% boxplot.stats(hw1_fj_dnaout]
hw1_fj_dna_IncomeFilter

不使用科學計數法

options(scipen=3)

20200925 14:01 RHW1###

1. 請將數據hw1_a和hw1_b分別讀入R论皆，查看數據并指出各個變量的形式益楼，最小值，最大值纯丸，中值偏形，均值静袖，標準差觉鼻。

讀取hw1_a.csv文件

install.packages("openxlsx")
library(openxlsx)
hw1_a<-read.csv("D:/mikiStudyFiles/les/hw1_a.csv")
hw1_a

方法1：陳肯-查看hw1_a各變量類型，最小值队橙，最大值坠陈，中值，均值捐康，標準差

str(hw1_a)

直接查看hw1_a中所有字段的類型

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

直接匯總求出hw1_a中所有字段的max min median mean 1st Qu. 3rd Qu.

方法2：miki

typeof(hw1_aAge)
typeof(hw1_aYears_at_Address)
typeof(hw1_aID,na.rm=TRUE)
max(hw1_aYears_at_Employer,na.rm=TRUE)
max(hw1_aIncome,na.rm=TRUE)
min(hw1_aAge,na.rm=TRUE)
min(hw1_aIncome,na.rm=TRUE)
min(hw1_aID,na.rm=TRUE)
median(hw1_aYears_at_Employer,na.rm=TRUE)
median(hw1_aYears_at_Address,na.rm=TRUE)
mean(hw1_aAge,na.rm=TRUE)
mean(hw1_aYears_at_Address,na.rm=TRUE)
mean(hw1_aID,na.rm=TRUE)
sd(hw1_aYears_at_Employer,na.rm=TRUE)
sd(hw1_aIncome,na.rm=TRUE)

讀取hw1_b.csv文件

hw1_b<-read.csv("D:/mikiStudyFiles/les/hw1_b.csv")
hw1_b

查看hw1_b各變量類型仇矾，最小值，最大值解总，中值贮匕，均值，標準差

typeof(hw1_bCredit_Card_Debt)
typeof(hw1_bIs_Default)
max(hw1_bCredit_Card_Debt)
max(hw1_bIs_Default)
min(hw1_bCredit_Card_Debt)
min(hw1_bIs_Default)
median(hw1_bCredit_Card_Debt)
median(hw1_bIs_Default)
mean(hw1_bCredit_Card_Debt)
mean(hw1_bIs_Default)
sd(hw1_bCredit_Card_Debt)
sd(hw1_bIs_Default)

2. 結合上課我們所學的幾種數據join 的形式花枫，將兩個數據集進行合并刻盐。對于每種數據合并的方式，請說明key, 并且報告合并后的數據樣本總行數劳翰。

left_join(data,by="連接字段")敦锌，nrow() 表示行數

left_join

hw1_lf<-hw1_a%>%
left_join(hw1_b,by="ID")
hw1_lf

right_join

hw1_rj<-hw1_a%>%
right_join(hw1_b,by="ID")
hw1_rj
nrow(hw1_lf)
nrow(hw1_rj)

inner_join

hw1_ij<-hw1_a%>%
inner_join(hw1_b,by="ID")
hw1_ij
nrow(hw1_ij)

full_join

hw1_fj<-hw1_a%>%
full_join(hw1_b,by="ID")
hw1_fj
nrow(hw1_fj)

semi_join

hw1_sj<-hw1_a%>%
semi_join(hw1_b,by="ID")
hw1_sj
nrow(hw1_sj)

anti_join

hw1_aj<-hw1_a%>%
anti_join(hw1_b,by="ID")
hw1_aj
nrow(hw1_aj)

3. 請篩選出hw1_a 中收入大于4000的樣本，并將此樣本和hw1_b 中Is_Default=1的樣本合并佳簸，你可以使用inner join的方式乙墙。這一問中你可以用pipe的書寫形式。

比較的時候一定要用==，=表示賦值

hw1_filter<-hw1_a%>%
filter(hw1_aIs_Default==1),by="ID"
)
hw1_filter

4. 在第2問的基礎上, 請給出Income對Years_at_Employer的散點圖听想，你發(fā)現了哪些趨勢和現象?

ctr+shif+c多行注釋

# #畫圖 ggplot(data,aes(x,y))+

# # geom_point()

na.rm=TRUE 排除空值

ggplot(hw1_fj,aes(Years_at_Employer,Income))+
geom_point(na.rm = TRUE)+
scale_y_continuous(labels = scales::comma)

5. 在第4問的基礎上 按照Is_Default 增加一個維度腥刹，請展示兩變量在不同違約狀態(tài)的散點圖。請使用明暗程度作為區(qū)分方式

alpha表示明暗程度 as.factor()表示轉換成向量

ggplot(hw1_fj,aes(Years_at_Employer,Income,alpha=as.factor(Is_Default)))+
geom_point(na.rm = TRUE)+
scale_y_continuous(labels = scales::comma)

6. 對于第5問汉买，請使用形狀作為另外一種區(qū)分方式肛走。

shape表示形狀

ggplot(hw1_fj,aes(Years_at_Employer,Income,shape=as.factor(Is_Default)))+
geom_point(na.rm = TRUE)+
scale_y_continuous(labels = scales::comma)

7. 請找出各個列的缺失值，并刪除相應的行录别。請報告每一變量的缺失值個數朽色，以及所有缺失值總數。

sum()表示計數组题，is.na()判斷是否為空值

方法1：陳肯-篩選出hw1_fj第2列不為空的所有記錄

hw1_fj=filter(hw1_fj,!is.na(hw1_fj[2]))

方法2：miki

計算缺失值總數

sum(is.na(hw1_fj))

計算每個變量的缺失值數

sum(is.na(hw1_fjAge))
sum(is.na(hw1_fjYears_at_Address))
sum(is.na(hw1_fjAutomobile_Debt))
sum(is.na(hw1_fjIs_Default))

查看缺失值的行

hw1_fj[!complete.cases(hw1_lf),]

刪除缺失值的行

hw1_fj_dna<-hw1_lf[complete.cases(hw1_lf),]

刪除第25列有空值的行

manghe[complete.cases(manghe[,25]),]

8. 找出Income中的極端值并濾掉對應行的數據

quantile()表示取百分位比的值

quantile(hw1_a$Income,c(0.025,0.975))
hw1_a2=filter(hw1_a,Income>14168.81&Income<173030.92)

先做直方圖

ggplot(hw1_fj_dna,aes(Income,))+
geom_histogram(bins = 30)+
scale_x_continuous(labels = scales::comma)

將頻數在1-5的數據放大

ggplot(hw1_fj_dna,aes(Income))+
geom_histogram(bins = 30)+
coord_cartesian(ylim = c(1,5))+
scale_x_continuous(labels = scales::comma)

通過圖形篩選出Income中的極端值

(hw1_fj_dna%>%
filter(hw1_fj_dnaIncome

刪除極端值行

hw1_fj_dna_IncomeFilter<-hw1_fj_dna[!hw1_fj_dnaIncome > 100000))$Income,]
hw1_fj_dna_IncomeFilter
nrow(hw1_fj_dna_IncomeFilter)

箱線圖連續(xù)型變量異常值判斷##

ggplot(hw1_fj_dna,aes(,Income))+

geom_boxplot()

去除極端值的箱線圖

boxplot(hw1_fj_dna$Income, outline = F)

boxplot(hw1_fj_dna$Income)

極端值

boxplot.stats(hw1_fj_dnaout

去除極端值

hw1_fj_dna_IncomeFilter<-hw1_fj_dna[!hw1_fj_dnaIncome)$out]

hw1_fj_dna_IncomeFilter

9. 將Income對數化葫男，并畫出直方圖和density curve.

方法1：陳肯-密度直方圖簡便畫法

inc<-hw1_a$Income
lninc<-log(inc)
hist(lninc,prob=T)
lines(density(lninc),col="blue")

方法2：miki

對數化Income直方圖

ggplot(hw1_fj_dna_IncomeFilter,aes(Income,))+
geom_histogram(na.rm = TRUE,bins=30,color="gray")+
scale_x_log10(labels = scales::comma)

對數化Income密度曲線

ggplot(hw1_fj_dna_IncomeFilter,aes(Income,))+
geom_density(na.rm = TRUE,color="pink",size=1)+
scale_x_log10(labels = scales::comma)

二個圖合在一個圖里很奇怪密度曲線會很平緩不明顯，這部分還需要琢磨琢磨

10. 以Income作為因變量崔列，Years at Employer作為自變量梢褐，進行OLS回歸，寫出回歸的方程赵讯，并指出自變量系數是否在某一顯著性水平上顯著盈咳。同時，解釋你的結果(這一問你自己發(fā)揮可以找code解決)边翼。

OLS回歸方程lm(y~x,data) summary()計算統計量

y<-hw1_fj_dna_IncomeFilterYears_at_Employer
miki<-lm(y~x,hw1_fj_dna_IncomeFilter)
miki

計算描述性統計量

summary(miki)

###################################20201007 RHW2####################################

1.編寫函數get.root(a,b,c),求解一元二次方程ax^2+bx+c=0的實根;

創(chuàng)建函數function(),if(){}else if(){}else{}

rt<-function(a,b,c){
if(a==0&b==0&c==0){

a,b,c均為0,等式恒成立鱼响，有無窮解,無意義

root<-NA
}
else if(a==0&b==0){

a,b均為0，c不為0,無解

root<-NULL
}
else if(a==0){

a為0,b,c不為0组底，有一個根

root<- -c/b
} else{

a,b,c均不為0丈积，判斷delta

d<-b^2-4ac
if(d>=0){
root<-(-b+c(1,-1)sqrt(d))/(2a)
} else
(root<-NULL)
}
return(root)}

2.已知一元二次方程ax^2+bx+c=0的三個系數都是隨機變量，其中a服從[1,5]上的均勻分布债鸡，b服從正態(tài)分布N(3,10),c服從均值為1的指數分布江滨，請編寫函數get.prop(),計算該方程有實根的概率；

rnorm(n, mean=0, sd=1) 高斯（正態(tài)）分布

rexp(n, rate=1) 指數分布

rgamma(n, shape, scale=1) γ分布

rpois(n, lambda) Poisson分布

rweibull(n, shape, scale=1) Weibull分布

rcauchy(n, location=0, scale=1) Cauchy分布

rbeta(n, shape1, shape2) β分布

rt(n, df) t分布

rf(n, df1, df2) F分布

rchisq(n, df) χ 2 分布

rbinom(n, size, prob)二項分布

rgeom(n, prob)幾何分布

rhyper(nn, m, n, k) 超幾何分布

rlogis(n, location=0, scale=1) logistic分布

rlnorm(n, meanlog=0, sdlog=1)對數正態(tài)

rnbinom(n, size, prob)負二項分布

runif(n, min=0, max=1)均勻分布

rwilcox(nn, m, n), rsignrank(nn, n) Wilcoxon分布

方法1：陳肯-使用循環(huán)計算

get.prob<-function(n){
a=runif(n,min=1,max=5)
b=rnorm(n,mean=3,sd=sqrt(10))
c=rexp(n,rate=1)
k=0
for (i in 1:n) {
if(sign(b[i]b[i]-4a[i]*c[i])==1|0)
{k=k+1}
}
return(k/n)
}

get.prob(100000)

方法2

創(chuàng)建函數function(),sum()計數 注意：N(3,10) 表示期望是3厌均，方差是10唬滑，標準差是根號10，所以表示為rnorm(n,3,sqrt(10))

p<-function(n){
a<-runif(n,1,5)
b<-rnorm(n,3,sqrt(10))
c<-rexp(n,1)
d<-b^2-4ac
sum(d>=0)/n
}

問題2：給定數據棺弊，請完成以下任務晶密，請給出code 和輸出結果。

1.請讀入數據镊屎，使用軟件分別給出 price, marketshare惹挟，和brand的缺失值數量。請按照每一個brand, 將數據按照先marjetshare 后price 進行從高到低排序缝驳；

讀數read.csv(),sum(is.na())計數缺失值

導入csv數據连锯，并統計缺失值數量

library(openxlsx)
library(tidyverse)
hw2<-read.csv("D:/mikiStudyFiles/les/hw2.csv")
hw2
sum(is.na(hw2marketshare))
sum(is.na(hw2$brand))

group_by()分組,summarise()分組摘要,arrange(desc())降序

按照brand分組归苍，并按marketshare,price降序排序

hw2_arr<-hw2[complete.cases(hw2),]%>%
group_by(brand)%>%
arrange(desc(marketshare,price))%>%
summarise(brand,marketshare,price)
hw2_arr

2.請按照brand 的種類，對price和marketshare 求均值运怖。

hw2_brand_avg<-hw2[complete.cases(hw2),]%>%
group_by(brand)%>%
summarise(avg_price=mean(price,na.rm=TRUE),avg_marketshare=mean(marketshare,na.rm=TRUE))
hw2_brand_avg

3.請按照brand 的種類拼弃，對price和marketshare 畫散點圖。

方法1：陳肯-分幾個小塊

ggplot(data=hw2[complete.cases(hw2),])+
geom_point(mapping = aes(x=price,y=marketshare))+
facet_wrap(~brand,nrow=2)

方法2：miki -按brand分顏色

hw2_point<-ggplot(hw2[complete.cases(hw2),],aes(price,marketshare,color=as.factor(brand)))+
geom_point(na.rm=TRUE)
hw2_point

4.請按照價格的均值摇展，產生新的變量price_new, 低于均值為“低價格”吻氧，高于均值為“高價格”。 同樣對市場份額也是咏连，產生變量marketshare_new, 數值為“低市場份額”和“高市場份額”

新增變量mutate(),將某個變量分組ifelse(條件,true,false)

price_avg<-mean(hw2[complete.cases(hw2),]$price,na.rm=TRUE)
hw2_mutate<- hw2[complete.cases(hw2),] %>%
mutate(
price_new = ifelse(
price >price_avg,'高價格','低價格'
),
marketshare_new=ifelse(
price >price_avg,'高市場份額','低市場份額'
)
)
hw2_mutate

case_when的用法

hw2_mutate <- hw2 %>%

mutate(

price_new = case_when(

price >price_avg ~ '低',

TRUE ~ '高'

))

5.請估計模型盯孙，marketshare為Y，price為X.

一般線性模型lm()

y<-hw2price
miki<-lm(y~x,hw2)

coef()表示調用函數

coef(miki)

summary()表示計算描述性統計量

summary(miki)

6. 請畫出（5）的擬合直線

方法1：陳肯-geom_abline()

ggplot(data=hw2)+
geom_point(aes(x=price,y=marketshare))+
geom_abline(data= m1,col= "blue")

方法2： miki- modelr::add_predictions()計算預測值

library(modelr)
hw2_mod<-hw2%>%
add_predictions(miki)
ggplot(hw2,aes(x),na.rm=TRUE)+
geom_point(aes(y=y),na.rm=TRUE)+
geom_line(aes(y=pred),na.rm=TRUE,
hw2_mod,
color="pink",
size=1)

7.請隨機產生若干直線祟滴，驗證（5）的結果是最優(yōu)的

方法1：陳肯-使用最小二乘法驗證擬合振惰，擬合準則是使 與 的距離 的平方和最小

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))

方法2：mikii-使用原理 取N條線，求最優(yōu)解垄懂，計算點到直線的距離最小的線

添加斜線函數geom_abline(),隨機取N個斜率和截距

models<-tibble(
k=runif(1000,-1,1),
b=runif(1000,-1,1)
)
ggplot(hw2,aes(x,y))+
geom_abline(aes(intercept=b,slope=k),models,alpha=1/4)+
geom_point(na.rm=TRUE)

optim()求最優(yōu)解

a1<-runif(1000,-1,1)
a2<-runif(1000,-1,1)
model1 <- function(a, hw2) {
a[1] + x * a[2]
}
model1(c(7, 1.5), hw2)
measure_distance <- function(mod, hw2) {
diff <- y[!is.na(y)] - model1(mod,hw2)[!is.na(model1(mod,hw2))]
sqrt(mean(diff ^ 2))
}
hw2_dist <- function(a1, a2) {
measure_distance(c(a1, a2), hw2)
}
models <- models %>%
mutate(dist = purrr::map2_dbl(a1, a2, hw2_dist))
models
best<-optim(c(0, 0),measure_distance,hw2)
best$par

計算每一個點到直線的垂直距離

計算預測值

model1<-function(k,b,hw2){
b+x*k
}

計算均方根誤差

measure_distance <- function(k1,b1,hw2) {
diff <- y[!is.na(y)] - model1(k1,b1,hw2)[!is.na(model1(k1,b1,hw2))]
sqrt(mean(diff ^ 2))
}

計算點到模型間的距離

library(purrr)
hw2_dist<-function(k,b){
measure_distance(k,b,hw2)
}
models<-models%>%
mutate(dist=map2_dbl(k,b,hw2_dist))
models

畫圖骑晶，取距離最小的10條線,顏色最亮的為最優(yōu)解

ggplot(hw2,aes(x,y))+
geom_point(size=2,na.rm=TRUE)+
geom_abline(aes(intercept=b,slope=k,color = -dist),
data=filter(models,rank(dist)<=10),na.rm=TRUE)

8.請估計模型，marketshare為Y草慧，price和brand 為X.

多變量線性模型lm(y~x1+x2+...,data)

y<-hw2price
x2<-hw2$brand
miki2<-lm(y~x1+x2,hw2)
coef(miki2)
summary(miki2)

#####################統計學基于R 20201017 8:30################

繪制核密度圖###

load("D:/mikiStudyFiles/les/example/ch2/example2_2.Rdata")
par(mfrow=c(1,2),cex=0.8,mai=c(0.7,0.7,0.1,0.1))
d<-example2_2$銷售額
plot(density(d),col=1,xlab="銷售額",ylab="密度",main="")
rug(d,col="pink")
plot(density(d),col="gold",border="black")
polygon(density(d),col="gold",border="black")
rug(d,col="brown")

用lattice包繪制核密度圖###

load("D:/mikiStudyFiles/les/example/ch2/example2_3_1.Rdata")
library(lattice)
dp1<-densityplot(~射擊環(huán)數|運動員,data=example2_3_1,col="blue",cex=0.5,par.strip.text=list(cex=0.6),sub="(a)柵格圖")
dp2<-densityplot(~射擊環(huán)數,group=運動員,data=example2_3_1,auto.key=list(columns=1,x=0.01,y=0.95,cex=0.6),cex=0.5,sub="(b)比較圖")
plot(dp1,split=c(1,1,2,1))
plot(dp2,split=c(2,1,2,1),newpage=F)

lattice包繪制的6名運動員射擊成績的和密度圖###

不同收入等級的臉譜圖###

注意安裝包括號內需要加雙引號install.packages(""),加載包不需要library()

install.packages("aplpack")
library(aplpack)
load("D:/mikiStudyFiles/les/example/ch2/matrix2_5.Rdata")
faces(matrix2_5,nrow.plot=4,ncol.plot=2,face.type=0)

計算繪制洛倫茨曲線所需的各百分比數值###

install.packages("DescTools")
library(DescTools)
Lc(example2_10人數)
Lc(example2_10人數)
Lc(example2_10人數)

繪制洛倫茨曲線###

plot(Lc(example2_10人數),xlab="人數比例",ylab="收入比例",col=4,panel.first=grid(10,10,col="gray70"))

ggplot繪圖###

三個不同專業(yè)的50名學生和統計學的考試分數數據###

set.seed(12345)
n=100
x<-round(rnorm(n,80,6))
y<-round(2+0.4*x+rnorm(n,50,3))
a<-sample(c("管理學","會計學","金融學"),20,replace=T)
b<-sample(c("男","女"),20,replace=T)
d<-data.frame(專業(yè)=a,性別=b,數學=x,統計學=y)
head(d)

20201102RHW3第三次作業(yè)################

問題 1：

A 和 B 約定在某籃球場見面桶蛔。他倆都不太守時，出現時間服從均勻分布漫谷。他倆也都沒有耐心仔雷， 每個人都會只等對方十分鐘就會離開。已知 A 到籃球場的時間為下午 4 點到 5點之間抖剿。

(1) 如果 B 到達籃球場的時間也為下午 4 點到 5 點之間朽寞，模擬運行 50000 次，看看他們成功相遇的概率斩郎。

方法1：陳肯-是否見面建一個表 tibble()表示簡單數據框,abs()表示絕對值

n_Sim <- 50000
sim_meet <- tibble(
A = runif(n_Sim, min = 0, max = 60),
B = runif(n_Sim, min = 0, max = 60)
) %>%
mutate(result = ifelse(abs(A - B) <= 10,
"They meet", "They do not"))
p_meet <- sim_meet %>% count(result) %>%
arrange(n) %>%
mutate(percent = n / n_Sim)
p_meet

方法2：miki

library(tidyverse)
a<-runif(50000,0,1)
b<-runif(50000,0,1)
c<-a-b
p<-sum(c>=-1/6&c<=1/6)
p/50000

(2) 對上一問的 50000 次模擬，用不同顏色在一張圖中展示成功相遇與否喻频。

d<-c>=-1/6&c<=1/6
hw3_1<-tibble(
A到達時間=a,
B到達時間=b,
AB到達時間差=c,
AB是否相遇=d
)
hw3_point<-ggplot(hw3_1,aes(A到達時間,B到達時間,color=AB是否相遇))+
geom_point()
hw3_point

(3) B 應該如何選擇 4 點到 5 點之間的哪個時間段缩宜，來提升他們成功相遇的概率？ 用模擬展示你的理由

分成6段甥温，看哪段的概率高

a1<-runif(50000,0,1)
b1<-runif(50000,0,1/6)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000
a1<-runif(50000,0,1)
b1<-runif(50000,1/6,1/3)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000
a1<-runif(50000,0,1)
b1<-runif(50000,1/3,1/2)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000
a1<-runif(50000,0,1)
b1<-runif(50000,1/2,2/3)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000
a1<-runif(50000,0,1)
b1<-runif(50000,2/3,5/6)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000
a1<-runif(50000,0,1)
b1<-runif(50000,5/6,1)
c1<-a1-b1
p1<-sum(c1>=-1/6&c1<=1/6)
p1/50000

問題 2： 請使用 nycflights13 和 pipe 語法

（1） 從 flights 數據表中挑選出以下變量：(year, month, day, hour, origin, dep_delay, distance, carrier)锻煌，將生產的新表保存為 flight1互广。

select(data,字段名1,字段名2...)

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

（2） 從 weather 數據表中挑選出以下變量： (year, month, day, hour, origin, humid, wind_speed)漫萄，將生產的新表保存為 weather1。

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

（3）將 flight1 表和 weather1 表根據共同變量進行內連接专缠，隨機抽取 100000 行數據狰挡， 將生產的結果保存為 flight_weather捂龄。 (提示：sample_n()函數释涛，不用重復抽取)

無放回隨機抽樣 sample_n(數據表,抽取行數,replace=Flase)

flight_weather<-sample_n(flight1%>%
inner_join(weather1,by=c("year","month","day","hour","origin")),100000, replace = FALSE)
flight_weather

（4）從 flight_weather 表中對三個出發(fā)機場按照平均出發(fā)延誤時間排降序，并將結果保 留在 longest_delay 表中倦沧。把結果展示出來唇撬。

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

（5）根據出發(fā)地（origin） 在同一個圖中畫出風速 wind_speed（x 軸）和出發(fā)延誤時間 dep_delay（y 軸） 的平滑曲線圖。

geom_smooth()表示平滑曲線

hw3_smooth<-ggplot(flight_weather,aes(wind_speed,dep_delay,color=origin))+
geom_smooth(se = FALSE,na.rm=TRUE)

se=False 去除陰影

hw3_smooth

（6）根據不同出發(fā)地（origin）在平行的 3 個圖中畫出風速 wind_speed（x 軸）和出發(fā) 延誤時間 dep_delay（y 軸）的散點圖展融。

將多個模型的可視化結果放在一張圖中 facet_wrap(~變量)

hw3_2_point<-ggplot(flight_weather,aes(wind_speed,dep_delay,color=origin))+
geom_point(na.rm=TRUE)+
facet_wrap(~origin)
hw3_2_point

（7）根據 flight_weather 表窖认，畫出每個月航班數的直方分布圖，x 軸為月份告希，y 軸是每個 月份航班數所占的比例扑浸。

geom_bar()表示條形圖,..prop..表示分類變量中每個類別占總量的比，group=1就是將這些類別當作一組的這樣一個整體去分別計算各個類別的占比

hw3_bar<-ggplot(flight_weather)+
geom_bar(aes(x=month,y= ..prop.., group = 1))
hw3_bar

（8）根據 flight_weather 表燕偶，畫出每個月航班距離的 boxplot 圖首装，x軸為月份，y 軸為 航行距離, 根據的航行距離的中位數從低到高對 x 軸的月份進行重新排序杭跪。

箱線圖

hw3_box<-ggplot(flight_weather,aes(x=month,y=distance,group = month))+
geom_boxplot()
hw3_box

根據中位數排序

reorder()重新排序

hw3_box_reorder<-ggplot(flight_weather)+
geom_boxplot(
aes(x=reorder(month,distance,FUN = median),y=distance
))
hw3_box_reorder

問題三：

（1)定義每個區(qū)域變種的熵為一定模式仙逻，熵是衡量事物分布均勻性的指標，寫出一個名為H的函數表示它

先建立矩陣

library(purrr)
a1<-c(0.2,0.2,0.2,0.2,0.2)
a2<-c(0.8,0.1,0.05,0.025,0.025)
a3<-c(0.05,0.15,0.7,0.05,0.05)
hw3_3<-rbind(a1,a2,a3)
hw3_3

給矩陣的字段命名

colnames(hw3_3) <- c("v1","v2","v3","v4","v5")
hw3_3

建立函數function(x,y,..,data)

H<-function(p,hw3_3){
-sum(p*log(p))
}
H

（2）定義區(qū)域與區(qū)間的距離為一個函數涧尿，寫一個名為KLD的函數表示它

KLD<-function(p,q,hw3_3){
sum(p*(log(p)-log(q)))
}
KLD
KLD(a1,a2)
KLD(a1,a3)
KLD(a2,a3)

（3） 應用purr：：map計算每個區(qū)域內變種分布的熵系奉，觀察并解釋各區(qū)該熵的大小順序

map(x,函數), list()函數合并多個列表成一個列表

hw3_map<-map(list(a1,a2,a3),H)

（4）計算每個區(qū)域與其他區(qū)域的KL距離，將結果保存為矩陣形式姑廉。該矩陣的可行可以代表每個區(qū)域到其他區(qū)域的變種分布的距離缺亮。

方法1：陳肯 先建立一個空矩陣，循環(huán)計算桥言，將計算得出的KLD值分別放入每個矩陣對應位置

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

方法2：miki-rbind(a,b)按行合并矩陣萌踱，cbind(a,b)按列合并矩陣

k1<-KLD(a1,a2)
k2<-KLD(a1,a3)
k3<-KLD(a2,a3)
hw3_KL<-rbind(k1,k2,k3)
hw3_KL

20201213 15:00統計學小組作業(yè)-盲盒######

詞云

library(wordcloud)
wordcloud<-read.csv("D:/mikiStudyFiles/les_business analyst/wordcloud.csv")
colors=c('red','blue','green','orange','MediumPurple','pink','DarkGray','blue','Gold')
wordcloud(wordcloudamount,scale=c(5,0.3),min.freq=-Inf,max.words=100,colors=colors,random.order=F,random.color=F,ordered.colors=F)

是否購買男女占比

library(tidyverse)
manghe<-read.csv("D:/mikiStudyFiles/les_business analyst/manghe.csv")
manghe<- manghe %>%
mutate(
是否購買 = ifelse(
是否購買過 == 0,'購買','未購買'
)
)

manghe[complete.cases(manghe[,25]),]
ggplot(manghe[complete.cases(manghe[,25]),]) +
geom_bar(aes(x = 性別, fill = 是否購買),na.rm=TRUE,position="dodge")

mh_new %>%
group_by(年齡,是否購買) %>% # calculate the counts
summarize(counts = n()) %>%
arrange(-counts) # sort by counts

是否購買年齡占比

ggplot(manghe[complete.cases(manghe[,25]),]) +
geom_bar(aes(x = 年齡, fill = 是否購買),na.rm=TRUE,position="dodge")

是否購買職業(yè)占比

ggplot(manghe[complete.cases(manghe[,25]),]) +
geom_bar(aes(x = 職業(yè), fill = 是否購買),na.rm=TRUE,position="dodge")

是否購買職業(yè)地域占比

ggplot(manghe[complete.cases(manghe[,25]),]) +
geom_bar(aes(x = 省, fill = 是否購買),na.rm=TRUE,position="dodge")

購買過盲用戶評價月消費

ggplot(manghe%>%filter(manghe是否購買過==1)) +
geom_bar(aes(x = 平均月消費),na.rm=TRUE)

抽到隱藏款概率

manghe<- manghe %>%
mutate(
是否抽到隱藏款 = ifelse(
是否抽到過隱藏款== 0,'未抽到過隱藏款','抽到過隱藏款'
)
)

ggplot(manghe[complete.cases(manghe[,26]),]) +
geom_bar(aes(x = 是否抽到隱藏款,y=..prop..,group=1),na.rm=TRUE)

上屆大神R：http://www.reibang.com/p/468862281e19