時系列データをあまり扱わないのでまじめに考えてなかったんですが、Rで単位根検定をする場合、拡張ディッキー＝フラー検定（augmented Dickey–Fuller test）を実施してくれるadf.test()という関数があります。
しかしこの関数は、
- 考慮するラグの次数を指定しなかった場合、サンプルサイズ（時系列データの長さ）を基準にして自動選択している
- 定数項もトレンド項も含むパターンしかやってくれない
という制限があり、後者についてはたいていの場合含めといたほうがいいらしいので問題ないとして、前者はよくわからない基準であって*1、本来はAICなどに基づいて選択されたラグ次数を用いたモデルで検定したほうがいいと思われます。実際ラグ次数によって結果がけっこう変わることもあるので、慎重に選んだほうがいいかと思います。
ということで、以下の記事で説明されているように、CADFtestパッケージのCADFtest()関数をつかっとけば簡単でいいんじゃないでしょうか。
上の記事だと色々やってて分かりづらいかもしれないので、パッケージの仕様書を見たほうがいいと思いますが、
これを読むと、共変量付きのADF検定ってのもできるらしいです。てうかそもそも，頭の"C"がCovariateのCで、それ用のパッケージみたいですね。
で、一番単純な使い方としては、
#dに時系列データがベクトルで入ってるとする CADFtest(d, type="trend", # トレンド項も定数項もあり max.lag.y=10, # ラグの最大次数を自分で適当に指定 criterion='AIC' # ラグ次数はAIC基準で選ぶ )
でいいんじゃないでしょうか（共変量を考えるときは、dを入れている最初の引数にモデル式を記述するとのこと）。
typeのところを"trend"にすると定数項もトレンド項もあり、"drift"にするとトレンド項なし、"none"にするとどっちもなしになるようです。よほどの理論的理由がない限り、両方入れておくのが保守的らしいです。
単純に結果だけ知りたい場合は、上記をprintして出てくるp-valueが0.05を下回っていれば、5%水準で「単位根あり」の帰無仮説を棄却できることになります。0.05より大きければ、単位根過程の疑いありです。
なお、正確に理解してないのですが、CADFtestでラグ0が選択されることがあり、そういう場合は、ADF検定ではなくDF検定に切り替えたほうがよく、adf.test()でk=0としてやればいいようです。（リンク）
詳しい使用例が以下のドキュメントにのってました。
Covariate Augmented Dickey-Fuller Tests with R
間違ってたらごめんなさい。
