It’s just playing around with quantile regression on a pooled dataset.

rm(list = ls())
setwd("/Users/vancam/Documents/WAIKATO-Thesis/VESuse/CVan/TFP_OwnGO")

load("/Users/vancam/Documents/WAIKATO-Thesis/VESuse/CVan/Final_forpaper/export_ublfinal.rda")
library(tidyverse)
## Loading tidyverse: ggplot2
## Loading tidyverse: tibble
## Loading tidyverse: tidyr
## Loading tidyverse: readr
## Loading tidyverse: purrr
## Loading tidyverse: dplyr
## Conflicts with tidy packages ----------------------------------------------
## filter(): dplyr, stats
## lag():    dplyr, stats
export <- filter(export_ublfinal, exportvol >0)

export <- mutate(export, loglabor = log(labor))
export <- mutate(export, size = labor)
export <- mutate(export, sizesq = labor*labor)
export <- mutate(export, logsize = log(size))
attach(export)
## The following object is masked _by_ .GlobalEnv:
## 
##     export
Y <- cbind(exportvol)
X <- cbind(size, scale, factor(ownership), Hlinkage, industrialzone)
library(quantreg)
## Loading required package: SparseM
## 
## Attaching package: 'SparseM'
## The following object is masked from 'package:base':
## 
##     backsolve
ols <- lm(Y~X, data = export)
summary(ols)
## 
## Call:
## lm(formula = Y ~ X, data = export)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -102717   -2104    -494    1098  394890 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      6567.1362   363.0450  18.089  < 2e-16 ***
## Xsize              11.0497     0.1427  77.412  < 2e-16 ***
## Xscale            243.4551     4.2749  56.950  < 2e-16 ***
## X               -4112.6279   309.4852 -13.289  < 2e-16 ***
## XHlinkage       -6248.6533   317.6572 -19.671  < 2e-16 ***
## Xindustrialzone  1253.8619   178.6808   7.017 2.34e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 10040 on 19161 degrees of freedom
## Multiple R-squared:  0.4292, Adjusted R-squared:  0.429 
## F-statistic:  2881 on 5 and 19161 DF,  p-value: < 2.2e-16
#Quantiles
quan10 <- rq(Y~X, data = export, tau = 0.1)
quan20 <- rq(Y~X, data = export, tau = 0.2)
quan30 <- rq(Y~X, data = export, tau = 0.3)
quan40 <- rq(Y~X, data = export, tau = 0.4)
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
quan50 <- rq(Y~X, data = export, tau = 0.5)
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
quan60 <- rq(Y~X, data = export, tau = 0.6)
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
quan70 <- rq(Y~X, data = export, tau = 0.7)
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
quan80 <- rq(Y~X, data = export, tau = 0.8)
quan90 <- rq(Y~X, data = export, tau = 0.9)
summary(quan10)
## Warning in summary.rq(quan10): 459 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.1, data = export)
## 
## tau: [1] 0.1
## 
## Coefficients:
##                 Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      15.24960   4.72373    3.22830   0.00125
## Xsize             0.16916   0.02985    5.66738   0.00000
## Xscale            1.06311   0.44298    2.39988   0.01641
## X                -3.58615   4.92986   -0.72743   0.46697
## XHlinkage       -10.21246   1.91136   -5.34302   0.00000
## Xindustrialzone   7.12891   3.45440    2.06372   0.03906
summary(quan20)
## Warning in summary.rq(quan20): 1146 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.2, data = export)
## 
## tau: [1] 0.2
## 
## Coefficients:
##                 Value     Std. Error t value   Pr(>|t|) 
## (Intercept)      78.06838   6.18869   12.61468   0.00000
## Xsize             0.97188   0.09234   10.52525   0.00000
## Xscale            3.10448   5.54116    0.56026   0.57531
## X               -60.00353   7.59780   -7.89749   0.00000
## XHlinkage       -38.76316   4.18777   -9.25627   0.00000
## Xindustrialzone  13.09452   5.48949    2.38538   0.01707
summary(quan30)
## Warning in summary.rq(quan30): 1402 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.3, data = export)
## 
## tau: [1] 0.3
## 
## Coefficients:
##                 Value      Std. Error t value    Pr(>|t|)  
## (Intercept)      211.55002   23.24330    9.10155    0.00000
## Xsize              2.11458    0.17570   12.03486    0.00000
## Xscale            25.18240   13.86862    1.81578    0.06942
## X               -173.05782   25.48325   -6.79104    0.00000
## XHlinkage       -110.87141    9.08251  -12.20713    0.00000
## Xindustrialzone    9.45440    2.65725    3.55796    0.00037
summary(quan40)
## Warning in summary.rq(quan40): 986 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.4, data = export)
## 
## tau: [1] 0.4
## 
## Coefficients:
##                 Value      Std. Error t value    Pr(>|t|)  
## (Intercept)      400.76055   33.85094   11.83898    0.00000
## Xsize              3.32785    0.27105   12.27753    0.00000
## Xscale            82.63497   31.40032    2.63166    0.00850
## X               -333.55664   33.05641  -10.09053    0.00000
## XHlinkage       -205.72697   15.88470  -12.95127    0.00000
## Xindustrialzone   10.29999   13.81507    0.74556    0.45594
summary(quan50)
## Warning in summary.rq(quan50): 815 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.5, data = export)
## 
## tau: [1] 0.5
## 
## Coefficients:
##                 Value      Std. Error t value    Pr(>|t|)  
## (Intercept)      474.90011   51.86717    9.15608    0.00000
## Xsize              4.34419    0.33993   12.77957    0.00000
## Xscale           231.47953   35.34815    6.54856    0.00000
## X               -371.36910   52.19331   -7.11526    0.00000
## XHlinkage       -309.53603   13.54900  -22.84567    0.00000
## Xindustrialzone   -6.43663    5.61223   -1.14689    0.25144
summary(quan60)
## Warning in summary.rq(quan60): 729 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.6, data = export)
## 
## tau: [1] 0.6
## 
## Coefficients:
##                 Value      Std. Error t value    Pr(>|t|)  
## (Intercept)      553.68819   82.60589    6.70277    0.00000
## Xsize              6.79834    0.42361   16.04876    0.00000
## Xscale           284.80673   32.20164    8.84448    0.00000
## X               -384.78890   83.04574   -4.63346    0.00000
## XHlinkage       -488.90205   24.24131  -20.16813    0.00000
## Xindustrialzone   13.52652   26.99804    0.50102    0.61636
summary(quan70)
## Warning in summary.rq(quan70): 377 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.7, data = export)
## 
## tau: [1] 0.7
## 
## Coefficients:
##                 Value      Std. Error t value    Pr(>|t|)  
## (Intercept)      771.85644   85.44348    9.03353    0.00000
## Xsize              9.64132    0.53984   17.85946    0.00000
## Xscale           384.62417   51.91733    7.40840    0.00000
## X               -472.17945   80.27040   -5.88236    0.00000
## XHlinkage       -754.30673   48.25191  -15.63268    0.00000
## Xindustrialzone   59.36857   45.24502    1.31216    0.18948
summary(quan80)
## Warning in summary.rq(quan80): 290 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.8, data = export)
## 
## tau: [1] 0.8
## 
## Coefficients:
##                 Value       Std. Error  t value     Pr(>|t|)   
## (Intercept)      1127.12436   164.83562     6.83787     0.00000
## Xsize              14.11493     0.62444    22.60432     0.00000
## Xscale            529.31077    70.39604     7.51904     0.00000
## X                -530.19854   160.72485    -3.29880     0.00097
## XHlinkage       -1255.08062    93.06621   -13.48589     0.00000
## Xindustrialzone   168.17722    94.16582     1.78597     0.07412
summary(quan90)
## Warning in summary.rq(quan90): 207 non-positive fis
## 
## Call: rq(formula = Y ~ X, tau = 0.9, data = export)
## 
## tau: [1] 0.9
## 
## Coefficients:
##                 Value       Std. Error  t value     Pr(>|t|)   
## (Intercept)      2396.02053   249.58022     9.60020     0.00000
## Xsize              19.86599     1.32495    14.99376     0.00000
## Xscale            853.12093   208.03324     4.10089     0.00004
## X                -786.65801   224.64980    -3.50171     0.00046
## XHlinkage       -2567.17855   166.23968   -15.44263     0.00000
## Xindustrialzone   391.77853    67.03016     5.84481     0.00000
 library(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2015). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2. http://CRAN.R-project.org/package=stargazer
stargazer(ols,quan10, quan20, quan30, quan40, quan50, quan60, quan70, quan80, quan90, title = "Export and firms' size - Quantile regression", type = "text")
## 
## Export and firms' size - Quantile regression
## ==============================================================================================================================================================
##                                                                                Dependent variable:                                                            
##                     ------------------------------------------------------------------------------------------------------------------------------------------
##                                                                                         Y                                                                     
##                                 OLS                                                                quantile                                                   
##                                                                                                   regression                                                  
##                                 (1)                 (2)        (3)         (4)         (5)         (6)         (7)         (8)          (9)          (10)     
## --------------------------------------------------------------------------------------------------------------------------------------------------------------
## Xsize                        11.050***            0.169***   0.972***   2.115***    3.328***    4.344***    6.798***    9.641***     14.115***     19.866***  
##                               (0.143)             (0.030)    (0.092)     (0.176)     (0.271)     (0.340)     (0.424)     (0.540)      (0.624)       (1.325)   
##                                                                                                                                                               
## Xscale                       243.455***           1.063**     3.104      25.182*    82.635***  231.480***  284.807***  384.624***   529.311***    853.121***  
##                               (4.275)             (0.443)    (5.541)    (13.869)    (31.400)    (35.348)    (32.202)    (51.917)     (70.396)      (208.033)  
##                                                                                                                                                               
## X                          -4,112.628***           -3.586   -60.004*** -173.058*** -333.557*** -371.369*** -384.789*** -472.179***  -530.199***   -786.658*** 
##                              (309.485)            (4.930)    (7.598)    (25.483)    (33.056)    (52.193)    (83.046)    (80.270)     (160.725)     (224.650)  
##                                                                                                                                                               
## XHlinkage                  -6,248.653***         -10.212*** -38.763*** -110.871*** -205.727*** -309.536*** -488.902*** -754.307*** -1,255.081*** -2,567.179***
##                              (317.657)            (1.911)    (4.188)     (9.083)    (15.885)    (13.549)    (24.241)    (48.252)     (93.066)      (166.240)  
##                                                                                                                                                               
## Xindustrialzone             1,253.862***          7.129**    13.095**   9.454***     10.300      -6.437      13.527      59.369      168.177*     391.779***  
##                              (178.681)            (3.454)    (5.489)     (2.657)    (13.815)     (5.612)    (26.998)    (45.245)     (94.166)      (67.030)   
##                                                                                                                                                               
## Constant                    6,567.136***         15.250***  78.068***  211.550***  400.761***  474.900***  553.688***  771.856***  1,127.124***  2,396.021*** 
##                              (363.045)            (4.724)    (6.189)    (23.243)    (33.851)    (51.867)    (82.606)    (85.443)     (164.836)     (249.580)  
##                                                                                                                                                               
## --------------------------------------------------------------------------------------------------------------------------------------------------------------
## Observations                   19,167              19,167     19,167     19,167      19,167      19,167      19,167      19,167       19,167        19,167    
## R2                             0.429                                                                                                                          
## Adjusted R2                    0.429                                                                                                                          
## Residual Std. Error   10,035.630 (df = 19161)                                                                                                                 
## F Statistic         2,881.167*** (df = 5; 19161)                                                                                                              
## ==============================================================================================================================================================
## Note:                                                                                                                              *p<0.1; **p<0.05; ***p<0.01
## Test for the difference in the intercepts between quantiles
anova(quan10,quan20, joint = FALSE)
## Warning in summary.rq(x, se = se, covariance = TRUE): 459 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 1146 non-positive fis
## Quantile Regression Analysis of Deviance Table
## 
## Model: Y ~ X
## Tests of Equality of Distinct Slopes: tau in {  0.1 0.2  }
## 
##                 Df Resid Df  F value    Pr(>F)    
## Xsize            1    38333 105.8018 < 2.2e-16 ***
## Xscale           1    38333   0.1507    0.6979    
## X                1    38333  68.8906 < 2.2e-16 ***
## XHlinkage        1    38333  55.8305  8.06e-14 ***
## Xindustrialzone  1    38333   1.7608    0.1845    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(quan10,quan30, joint = FALSE)
## Warning in summary.rq(x, se = se, covariance = TRUE): 459 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 1402 non-positive fis
## Quantile Regression Analysis of Deviance Table
## 
## Model: Y ~ X
## Tests of Equality of Distinct Slopes: tau in {  0.1 0.3  }
## 
##                 Df Resid Df  F value    Pr(>F)    
## Xsize            1    38333 138.5613 < 2.2e-16 ***
## Xscale           1    38333   3.1215   0.07727 .  
## X                1    38333  51.4197 7.594e-13 ***
## XHlinkage        1    38333 141.1708 < 2.2e-16 ***
## Xindustrialzone  1    38333   0.4265   0.51373    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(quan10,quan20,quan30,quan40,quan50,quan60,quan70,quan80,quan90, joint = FALSE)
## Warning in summary.rq(x, se = se, covariance = TRUE): 459 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 1146 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 1402 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 986 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 815 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 729 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 377 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 290 non-positive fis
## Warning in summary.rq(x, se = se, covariance = TRUE): 207 non-positive fis
## Quantile Regression Analysis of Deviance Table
## 
## Model: Y ~ X
## Tests of Equality of Distinct Slopes: tau in {  0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9  }
## 
##                 Df Resid Df F value    Pr(>F)    
## Xsize            8   172495 68.4291 < 2.2e-16 ***
## Xscale           8   172495 14.2249 < 2.2e-16 ***
## X                8   172495 22.8420 < 2.2e-16 ***
## XHlinkage        8   172495 88.0873 < 2.2e-16 ***
## Xindustrialzone  8   172495  7.8102 1.523e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Ploting the data
kevin <- rq(Y~X, tau = seq(0.05, 0.95, by = 0.05), data = export)
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique

## Warning in rq.fit.br(x, y, tau = tau, ...): Solution may be nonunique
kevin.plot <- summary(kevin)
## Warning in summary.rq(xi, U = U, ...): 459 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 1190 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 1146 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 1240 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 1402 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 989 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 986 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 749 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 815 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 987 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 729 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 573 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 377 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 282 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 290 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 164 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 207 non-positive fis
## Warning in summary.rq(xi, U = U, ...): 2 non-positive fis
plot(kevin.plot)

phanvi<- rq(exportvol ~ size, data = export, tau=c(0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9))

van <- data.frame(t(coef(phanvi)))
colnames(van) <- c("intercept", "slope")
van$quantiles <- rownames(van)
van
##           intercept      slope quantiles
## tau= 0.1   7.077566  0.1982763  tau= 0.1
## tau= 0.2   2.019794  0.9982005  tau= 0.2
## tau= 0.3  -6.844961  2.2635659  tau= 0.3
## tau= 0.4 -14.421348  3.9044944  tau= 0.4
## tau= 0.5 -21.714286  6.1428571  tau= 0.5
## tau= 0.6 -34.182082  9.3986217  tau= 0.6
## tau= 0.7 -38.663805 13.8784514  tau= 0.7
## tau= 0.8  34.365471 19.8049327  tau= 0.8
## tau= 0.9 373.484211 33.2315789  tau= 0.9
p<- qplot(x=size, y=exportvol, data=export)
p+ geom_abline(aes(intercept=intercept, slope=slope,colour=quantiles), data=van)